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We present a new analytic tool for the study of cosmic structure formation, a 
double distribution of the number density of dark matter halos with respect to both 
halo mass and local over- (or under-) density. The double distribution provides a 
statistical treatment of the properties of matter surrounding collapsed objects, and 
can be used to provide analytical insight into environmental effects on hierarchical 
structure formation. 

We apply this new tool to the case of cosmic accretion shocks. We investigate 
and quantify the effect of environmental factors on the statistical properties of these 
shocks. For this purpose, we explore two different models. The first "control" model 
uses a Press-Schechter mass function to describe the population of collapsed struc- 
tures, and assumes that all objects accrete gas of the same density and temperature. 
The second model treats the accreted material as a multi-temperature, multi-density 
medium with densities and temperatures derived from the double distribution. We 
find that the shock environment significantly alters the physical impact of cosmic 
accretion shocks on the intergalactic medium, as well as the cosmic history of their 
properties. 
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We present a new analytic tool for the study of cosmic structure formation, a 
double distribution function of dark matter halos, with respect to both object mass 
and local over- (or under-) density. The double distribution provides a statistical 
treatment of the properties of matter surrounding collapsed objects, and can be used 
to study environmental effects on hierarchical structure formation. 

We apply this new tool to the case of cosmic accretion shocks. We investigate 
and quantify the effect of environmental factors on the statistical properties of these 
shocks. For this purpose, we explore two different models. The first "control" model 
uses a Press-Schechter mass function to describe the population of collapsed struc- 
tures, and assumes that all objects accrete gas of the same density and temperature. 
The second model treats the accreted material as a multi-temperature, multi-density 
medium with densities and temperatures derived from the double distribution. We 
find that the shock environment significantly alters the physical impact of cosmic 
accretion shocks on the intergalactic medium, as well as the cosmic history of their 
properties. 



Chapter 1 
Introduction 



The large-scale structure of the universe has been investigated with increasing 
intensity for almost a century. Cosmologists have been using multi-wavelength obser- 
vations, analytical models, and massive numerical simulations to study the properties 
of galaxies, galaxy clusters, and superclusters, their formation, and their evolution. 

Structure in the universe is believed to be seeded by small, random density fluctua- 
tions in the early universe, which grow gravitationally in a process further complicated 
by the interplay between baryons, radiation, dark matter, and dark energy. In this 
context, the observationally relevant predictions of any theory of structure forma- 
tion are the statistical properties of the resulting, observable population of cosmic 
structures. 

1.1 The Double Distribution 

Cosmological distributions have long been used for the study of cosmic structure 
formation with great success as an analytical tool complementary to numerical sim- 
ulations. They have been used to constrain the cosmological parameters; interpret 
results of cosmological simulations; study regions of the parameter space which can- 
not be approached by simulations due to prohibitive computational cost; explore the 
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effects of various physical processes in an efficient if approximate way. The analyt- 
ical tool used most widely in cosmology is the mass function of dark matter halos 
(distribution of the number density of halos with respect to halo mass). Analytical 
descriptions of dark matter halos are usually based on the Press-Schechter formal- 
ism ([1; 2]) and its extensions (e.g., [3-6]). The Press-Schechter mass function has 
been shown to agree well with results of N-body simulations (e.g., [7-9]). More so- 
phisticated approaches taking into account deviations from spherical symmetry (e.g., 
[10-12]) have improved this agreement even further. 

To derive the Press-Schechter mass function, one begins with a density perturba- 
tion field still in its linear regime. Regions in space are then smoothed on successively 
smaller scales. The mass of a collapsed object is then taken to be the largest smooth- 
ing mass scale for which the average linear overdensity exceeds some threshold. In 
this way, matter in the universe is distributed among collapsed structures of different 
masses, which all share the same value of average overdensity (the threshold value). 
Information about the local environment of collapsed objects (whether they live in 
underdensities or overdensities) is thus erased. 

For this reason, and despite its wide applicability, the expression for the mass 
function cannot be used to address environment-related questions: Does the mass 
function of structures in superclusters differ from the mass function inside voids? 
Are structures of some particular mass more likely to reside inside underdense or 
overdense regions in space? How does such a preference evolve with redshift, and 
how sensitively does it depend on the cosmological parameters? How does the state 
of the material surrounding and accreted by a collapsed object depend on the mass 
of the object and the cosmic epoch? 

To address such questions, we seek a double distribution of the number density 
of structures with respect to mass but also to local overdensity (or underdensity) . 
In order to extract information about the surroundings of collapsed structures, we 
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use the same random walk formalism which rigorously yields the Press-Schechter 
mass function ([2; 3]). Integration of this distribution over density contrast should 
return the Press-Schechter mass function so that the successes of the Press-Schechter 
formalism be retained. 

Two complications arise in the effort to expand the Press-Schechter mass function 
to incorporate a description of the local overdensity. First, the concept of the "local 
environment" is somewhat vague and needs to be defined in a more rigorous way. The 
size of the "local environment" cannot be the same for all structures. If this was the 
case, very small structures would represent only a tiny fraction of the "environment" , 
while very large structures could even exceed the size of the "environment" , which 
would be an unphysical situation. This problem is not exclusive to analytical tools, 
but also needs to be addressed when analyzing the results of numerical simulations. 
Second, the Press-Schechter treatment of the density field uses linear theory, and 
ways of converting this information to a more physical non-linear result need to be 
determined. 

We address the first problem by introducing a clustering scale parameter, j3, which 
allows us to define the size of the "environment" of each structure as a function of its 
mass. We address the second concern by calculating conversion relations between the 
linear-theory overdensities (or underdensities) and those predicted by the spherical- 
evolution model. 

Environment-related questions in cosmological structure formation have also been 
addressed using analytical models for the clustering properties of dark halos which 
evaluate quantities such as the cross-correlation function between dark halos and 
matter, and the biasing factor (e.g. [13-19], also see review by [20] and references 
therein). These analyses are based on the same random walk formalism which we use 
here to derive our double distribution (also see [21] for fitting formulae from N-body 
simulation results). However, the information content of the double distribution, 
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which treats the "environment" in a mass-dependent fashion, is complementary to 
that of correlation functions, which describe the clustering properties of the dark 
halo population at some fixed spatial scale. The double distribution is ideally fitted 
for population studies of cosmological objects. If the properties of a single object 
can be parametrized as a function of its mass and its environment, then the double 
distribution can be used to predict the statistical properties of such objects, as well 
as their evolution with time, for any cosmological model. 

We present such an application of the double distribution for the case of cosmic 
accretion shocks. The properties of cosmic accretion shocks depend sensitively on 
their environment, and for this reason the double distribution is the tool of choice to 
study this population of cosmic structures. 

1.2 Cosmic Accretion Shocks 

The formation of shocks in the baryonic component of matter in the universe is 
an inevitable and integral part of the process of cosmological structure formation. 
The longstanding question of the dark baryons 1 [22-25] is likely on its way to res- 
olution with the first evidence [26-29] for a large component of diffuse, low-density 
intergalactic gas distributed in the filaments that comprise the "cosmic web," and 
detected in the X-ray forest [30-32] via absorption lines from highly-ionized metals. 
This warm-hot intergalactic medium is thought to arise in structure formation shocks 
[30; 33-36]. Moreover, structure formation shocks heat the intergalactic medium and 
are likely to act as acceleration sites for nonthermal particles, the "structure forma- 
tion cosmic rays" (e.g. [37-43] and references therein). Such a cosmic-ray population 
would have distinct 7-ray and light-element signatures, which are currently subjects 
1 & fraction of the baryonic content of the universe which, until recently, has eluded detection, 
presumably because its temperature (10 5 — 10 7 K) is too high to produce Lya absorption while its 
density is too to be detected in emission with currently available equipment 
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of intense investigation [40; 44-58]. 

Cosmic shocks occur during different facets of structure formation and in a variety 
of environments, hence there are two possible ways to categorize them: according to 
the physical processes causing them, and according to the state of the medium in 
which they form. 

There are three principal processes associated with cosmic structure formation 
which result in the formation of large-scale shocks. 

1. Accretion of intergalactic gas by a collapsed, virialized structure. In this case, 
an accretion shock is formed at the interface between virialized and diffuse gas 
[59-61]. The shock is driven by the gravitational attraction exerted on the 
diffuse gas by the accretor. 

2. Merger of two collapsed structures. In this case, a merger shock is formed at 
the interface between the gas components of the merging objects [41; 61]. The 
shock is driven by the mutual gravitational attraction between the objects. 

3. Blast-like expansion of a void (an underdensity in its non-linear regime). The 
regions compressed between expanding voids form large-scale filaments, and 
filament shocks are formed at the interface between the expanding void and the 
compressed gas [62]. In this case, the shocks are driven by the expansion of the 
void rather than the gravity of the filaments. 

Shocks can also be divided according to the state of the gas passing through 
them, into external and internal shocks [63]. External shocks process pristine mate- 
rial, which has never been shocked before by any of the processes described above. 
External shocks are mostly filament shocks, since the process of formation of indi- 
vidual virialized structures (associated with the other two types of shocks) occurs 
principally within filaments, and therefore in most cases involves gas which has al- 
ready been processed at least by filament shocks. Internal shocks process gas which 
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has already been shock- heated in the past. All merger shocks, as well as many accre- 
tion shocks, are internal shocks. 

Since external shocks process colder material of lower sound speed, their Mach 
numbers are generally higher than those of internal shocks. However, because the gas 
passing through internal shocks has already been compressed, internal shocks process 
more mass and kinetic energy than external shocks. 

In the second part of this thesis, we present an analytical study of the population 
of accretion shocks, in a concordance Q m + Q\ = 1 universe. We use the double 
distribution of cosmic structures to describe the accreting structures as well as their 
environment at a given redshift. Other Press-Schechter extensions have been used to 
model different "families" of large-scale cosmic shocks by [41] and [64], who model 
merger shocks, and by [35], who describe large-scale shocks that may appear when 
overdense perturbations reach and exceed their turnaround point. Here, we present 
the first analytic model for accretion shocks, and we include, for the first time, a 
detailed treatment of the environment in which accreting structures reside. 

This thesis is organized as follows. In chapter 2 we derive the double distribution of 
cosmic structures, and we explore its information content by plotting the distribution 
itself as well as interesting derivative quantities for different cosmological models. In 
chapter 3 we apply the double distribution to study the effect of environment on the 
statistical properties of the population of cosmic accretion shocks around collapsed 
structures. We discuss our findings in chapter 4. 
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Chapter 2 

The Double Distribution Of 
Cosmic Structures 

2.1 Overview 

The double distribution of cosmic structures with respect to mass and local over- 
density is a cosmological statistical distribution which describes how the number 
density of collapsed and virialized dark matter objects is distributed among different 
masses and among different local density contrasts with respect to the cosmic mean 
density. In this chapter we derive the double distribution of cosmic structures and 
we explore its information content. 

The derivation of the double distribution is based on the "random walk" formal- 
ism, which was introduced for the derivation of cosmological mass functions from an 
early (linear) field of density fluctuations by [2] and by [3]. The basic simplifying 
assumptions behind building a mass function starting from a linear overdensity field 
are: 

1. All matter in the universe can be distributed among collapsed objects, and there 
is no diffuse matter in the universe. 
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2. At any given cosmic epoch, all structures can be viewed as though they have 
just virialized. 

3. All collapsed objects at a given cosmic epoch have the same mean matter den- 
sity 

4. All collapsed objects at a given cosmic epoch are spherically symmetric. 

5. The mean density of a certain region at a very early cosmic epoch (while all 
density perturbations are still growing linearly) contains adequate information 
to predict the mean density of the same region at much later times. 

The conceptual idea behind the process necessary to produce a mass function is the 
following. 

• We focus at a infinitesimal mass element centered on a fixed point in space. 

• Starting at an infinite mass scale and gradually proceeding to smaller mass 
scales, we smooth the linear density field and evaluate the mean overdensity 
at the scale under consideration (note that at an infinite mass scale the mean 
overdensity is always zero, since the smoothing identically returns the mean 
cosmic density). 

• We convert the smoothed linear overdensity to a nonlinear overdensity using 
the spherical evolution model. If the result is equal to the mean nonlinear 
overdensity of virialized objects, then the mass scale at which the smoothing 
took place is the mass of the collapsed object which hosts the infinitesimal mass 
under consideration. If not, we repeat the smoothing at a smaller mass scale, 
until the desired mean overdensity is reached. 

• We then repeat the process for all points in space, and thus calculate how much 
mass corresponds to every collapsed object mass interval. 
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The random walk formalism now recognizes that the smoothed overdensity executes 
a "random walk" with changing smoothing mass scale, and uses this fact to calculate 
the probability that any given point in space belongs to a collapsed object with mass 
in a specific interval. This probability can then be manipulated to give us the mass 
function. The mathematical details of this process are briefly reviewed in section 2.2. 

In our case, we would like to additionally retain information about the density of 
the region within which a collapsed structure is embedded. For this reason, we need 
to keep track of the mean density evaluated at a random walk "step" (a smoothing 
mass scale) preceding the final step which returns the virial density of the collapsed 
object. The "distance" in mass between the two "steps" of interest will determine 
how far from the collapsed object the "local environment" extends. In our derivation 
of the double distribution the ratio of the mass scale representing the environment 
of an object (including the object itself) over the mass of the object is defined as fl, 
the "clustering scale parameter" . It is a free parameter in our model, and, for every 
realization of our model, it is the same for all structures. In section 2.3 we present 
the details of the derivation of the double distribution through this extension of the 
random walk formalism. 

As we discussed above, a critical step in deriving both the Press-Schechter mass 
function as well as the double distribution of cosmic structures is the conversion 
of a density contrast evaluated using linear theory to a density contrast evaluated 
using an exact (if idealized) model. Since we have assumed that all structures are 
spherically symmetric, we will consider the "true" density contrast to be the one given 
by the spherical evolution model. In section 2.4 we derive such conversion relations 
between the linear and the spherical evolution models, by expanding the spherical 
model around zero overdensity and demanding that at early times the growth of 
density perturbations according to both models agree to first order. We derive exact 
conversion relations for both Einstein-de Sitter (flat, matter-only) and ACDM (flat, 
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matter + cosmological constant) universes. In addition, we present a useful and simple 
approximate conversion relation, which exhibits the correct asymptotic behavior at 
early as well as late times, and has an accuracy better than 2% throughout its domain. 

In section 2.5 we define several interesting integral moments of the double dis- 
tribution, describing the statistical properties of the population of collapsed cosmic 
structures and their distribution among different (overdense or underdense) environ- 
ments. Finally, in section 2.6 we explore the information content of the double distri- 
bution by plotting both the double distribution itself as well as several of its integral 
moments, for different redshifts and different cosmologies. Our principal finding is 
that for every cosmic epoch and every collapsed object mass, there is a pronounced 
peak in the double distribution, corresponding to a "most probable" overdensity (or 
underdensity) for the local environment of the structure. The location of this peak 

• moves towards higher overdensities with increasing redshift; 

• moves towards higher overdensities with increasing object mass; 

• depends only very mildly on the value of f3 for low j3 (while it eventually moves 
towards zero overdensity as (3 — * oo) 

• has a qualitatively similar behavior for both cosmological models we examined. 
At the present cosmic epoch, most structures are located inside underdensities. 

2.2 Random Walks and the Press-Schechter Mass 
Function 

The Press-Schechter mass function of collapsed structures is the comoving number 
density of virialized objects per differential mass interval, dn/dm, for every cosmic 
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epoch a 1 . A related quantity is the mass fraction, P(> m, a), which is the fraction of 
matter in the universe belonging to collapsed structures with mass > m. If P(> m,a) 
is known, then dn/dm can be calculated from 



where p m> o is the present-day matter density of the universe. 

P{> m, a) is in turn calculated by assigning, at each epoch a, every infinitesimal 
element dm in the universe to a collapsed structure of some mass m. A structure is 
considered "collapsed" if its mean overdensity 



exceeds a certain critical value, S c (a). In Eq. (2.2), p ma is the mean matter density 
of the universe at epoch a. The critical overdensity 5 c (a) is the mean overdensity 
predicted by the spherical evolution model for a structure virializing at epoch a. 
For each point in space, one calculates the mean local overdensity by smoothing the 
overdensity field 5(x,a) with a spherically symmetric filter function of varying mass 
scale, starting from m — > oo, where one averages over the whole universe and finds 
identically (5) = 0, and proceeding to successively smaller scales. When a mass scale 
is found for which the mean overdensity becomes equal to S c (a), it is taken to be the 
mass m of the parent object of the infinitesimal mass at the point under consideration. 
This way of assigning object masses circumvents the structure-in-structure problem, 
since the mass of the parent object is always the largest possible mass satisfying the 
criterion for collapse. All information on substructure within collapsed structures is 
thus erased from the resulting mass function. 

The way the average overdensity (5) changes when the smoothing mass scale is 
lr rhe proper mass function, i.e. the number of collapsed structures per unit proper volume 
per differential mass interval is simply related to the comoving mass function via dn/dm | pro per — 

d dll I dm | comoving 





'structure/ Pm,a 




Pm,a 
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varied resembles, under certain conditions, a ID random walk [2]. For all "parti- 
cles" (points in space in our case), the walk begins at the "spatial origin" ((5) = 0), 
at "time zero" (m — > oo). As "time progresses" (m decreases), each "particle" may 
move either to the "left" (negative (5)) or to the "right" (positive (5)). An "absorbing 
wall" exists at 5 c (a). If this "wall" is reached, the "particle" is "removed" from the 
walk (the point is assigned its parent object mass and removed from further consid- 
eration at smaller values of m). P(> m,a) is then the fraction of "particles" which 
have been "lost" by "time" m, and it can be calculated using random-walk theory. 

However, we must first ensure that simple random- walk theory is indeed appli- 
cable. First, each "step" of the "walk" should be completely independent from the 
previous step. This requires that the /c-modes producing an increase A (5) in the 
space-like variable not appear in any of the previous steps in (5). A smoothing win- 
dow function sharp in k— space, 

{1 k < kJm) 
(2.3) 
k > k c {m) 

(see [2] and [3] for more extended discussions on the consequences of such a choice) 
enforces this condition, since 

(<Wo = fw m (\x -x\) 5(x)d 3 x = Le^d 3 k (2.4) 

J Jk<k c (m) 

and 

A(S)g = / 6 k e*°- n <Pk (2.5) 

kc (jti) < fc< fcc (m — dm) 

which only involves new fc-modes corresponding to scales from m to m — dm. 

Second, there must be an equal probability for the system to "move" towards any 
one of the two available "directions" . A Gaussian overdensity field (which is the usual 
assumption for deriving analytic mass functions and which we adopt here) guarantees 
that this condition is satisfied. 
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Finally, the appropriate "time-like" variable (which should depend on m) needs 
to be selected, given that the "space-like" variable is (5). By direct analogy to the 
ID random walk theory result (x 2 ) = 2Dt, and from the definition of the variance of 
the overdensity field S(m), 

rk(m) 

S(m) = a 2 (m) = (\S(m,x)\ 2 ) oc / k 2 dk\S k \ 2 (2.6) 

Jk=0 

we can immediately identify Dt — > S(m)/2. 

Three further complications need to be addressed. First, our knowledge of 5k 
and subsequently S(m) is limited at late times. In the early universe, right after 
matter-radiation equality, (|<5fc| 2 ) can be simply described in terms of a power-law in 
k modified by a transfer function, (|5fc| 2 ) oc T 2 (k)k n . While all S are still in their linear 
regime, they simply grow by the linear growth factor (independent of k). However, 
at later times, when certain structures start departing from the linear regime, we 
cannot use our simple early-universe expressions for 8^. Second, (5) is limited to be 
> —1, which introduces a second, reflecting "wall" at a value of (5) = — 1, further 
complicating the random-walk calculations. Finally, the true overdensity field loses 
its Gaussianity as it evolves past the linear regime. 

To circumvent these problems, we define the linearly extrapolated overdensity field, 
S(x,a), as the overdensity field that would result if all structures continued to grow 
according to the linear theory until time a. Now S(x,a) is not limited to be > — 1, 
since it does not represent real overdensities. In addition, we can always calculate 
S(m) for 5(x,a), since 8k is modified from its simple early-universe expression only 
by the linear growth factor. Finally, the extrapolated field remains Gaussian at all 
times. 

The linearly extrapolated overdensity 5(x,a) and the associated variance S(m), 
are time- varying, but the time dependence is well-known (see chapter 3), and the 
same for both S and 5 2 2 . Thus, the time dependence drops out of ratios 8/y/S which 
2 as seen by Eq. (2.6) re-written for the extrapolated rather than the true overdensity field 
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appear in the mass function. For this reason we may, without loss of generality, choose 
any single epoch to evaluate these quantities, with the stipulation that 5 and S(m) 
must refer to the same epoch. Given that a(m) is often normalized to the present 
value of as, a convenient choice of epoch is the present. Then, Eq. (2.6) gives for 
S(m) 

8 /'i™* 1 'T*(k)k"+*dk' 
Thus we only consider S(x, ao) (the overdensity field linearly extrapolated to the 
present epoch), which we use instead of the true field S(x, a) in our random walk 
formalism 3 . To find the mass function at a particular cosmic epoch a, we calculate 
the location of the "absorbing wall", 5 c (a). If a structure is predicted to collapse at 
epoch a according to the spherical evolution model, then <5o jC (a) is the overdensity this 
same structure would have had if, instead of turning around and collapsing, it had 
continued its linear evolution until the present. This 5o,c(«) is then our "absorbing 
wall". 

We can now derive the mass fraction and mass function using random walk theory. 
If a particle executes a one-dimensional random walk with an absorbing boundary at 
a point x±, then its probability W(x, t) to be between x and x + dx at time t is [65] 

(2xi -x) 2 ' 



cxp 



X 2 



exp 



4Dt 



W(x, t, x\)dx = - — -. -dx , (2 



where x < x\. In our case, the probability that a point in space will be assigned 
an average extrapolated overdensity between 5 and S + dS when filtered at a scale m 



3 Physically, the substitution of the true field by the extrapolated field in the "random walk" 
corresponds to smoothing the extrapolated overdensity field, and then mapping the mean extrapo- 
lated overdensity value to a true overdensity value. That true overdensity value is then assumed to 
accurately represent the result of a smoothing of the true field, which implies 6 (J^/J — (^))- 
This would be exactly true only if 6(5) was linear in 6, which is not the case (see chapter 3). This 
assumption introduces an inaccuracy inherent to all calculations which employ it, including the 
Press-Schechter mass function as well as the double distribution. 
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corresponding to a variance of S(m) is 

exp 

W(S, S,S , c )dS = 



n 

'25 



exp 



(26q, c -S) 2 
25 



dS . 



(2.9) 



with 5 < 5o, c - The mass fraction P(> m, a) is then the fraction of points already 
"lost" from the walk when filtering at higher mass scales, which is one minus the 
fraction of points remaining in the walk, 



P(> m, a) 



P(>8 



0,c) 



1 - 



So.c 



W(5,S,8 ,c)dS~ 



erfc 



( <W(a) 



Then, 



dP(> m, a) 
dm 



1 5 0c (a) dS 



2n S(m) 3 / 2 dm 



exp 



2S{m) 



and the Press-Schechter mass function can be found using Eq. (2.1) 



dn 
dm 



m, a) 



2 p m ,o ^o,c(a) 
y/S(m 



it rrr 



dln^S 



dhim 



exp 



£q,c(q) S 
' 2S(m) 



(2.10) 



(2.11) 



(2.12) 



2.3 Derivation of the Double Distribution 

We now use the random walk formalism described in the previous section to derive 
the double distribution of the comoving number density of collapsed structures with 
respect to object mass m and local environment overdensity 5g, dn/(dmd5e). 

For the reasons described in the previous section, we replace the true overden- 
sity field, 8(x, a), with its linear extrapolation to the present time, S(x,a ). Thus, 
we derive the double distribution of comoving n with respect to object mass m and 
extrapolated local environment overdensity Sg, dn/(dmd8e). We then use linear the- 
ory and the spherical evolution model to establish a conversion relation 8(5, a) and 
calculate dn/(dmd5i) as 



dn 



dmd5. 



■(St, m, a)dmdSi 



dn 



dmdSn 
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Si(St, a),m,a 



dm^—dSz . 

OQp 



(2.13) 



First of all, we need to define the local environment extrapolated overdensity 
Si in a precise way. We would like 5i to be a measure of the density contrast of 
the medium in which a collapsed structure is embedded. Clearly, the value of 8g 
depends on how far from the structure itself its "environment" extends. We quantify 
this notion by introducing the clustering scale parameter, j3, which is defined in the 
following way: the "environment" of an object of mass m is a surrounding region in 
space which encompasses mass (3m (including the mass of the object). Hence, the 
local environment extrapolated overdensity 5e is the result of a filtering of 5(x, do) 
with a filter of scale (3m centered on the object. 

Formally, the above definition translates as follows. Consider the sharp in k— space 

filtering function W m {k) discussed previously (Eq. 2.3). The relation between the 

cutoff wavenumber and the filter mass, k c {m), is found by considering the form of the 

filter function in configuration space, 

sin [k c (m)r] - k c {m)r cos [k c (m)r] 
vv m [r) - , (^.14) 

and multiplying by p mi0 and integrating over all space, which yields [3], 

Um)= (^£]'\ (2.15) 



m 

For a collapsed structure at an epoch a which has mass m and is centered at a point 
Xo, we can write 

6(m,x ) = J W m (\x - x\) 6(x, a )d 3 x = S 0>c {a) , (2.16) 

since the mean extrapolated overdensity of the collapsed structure itself is always 
the critical value for collapse, So jC (a). For that same object, the local environment 
extrapolated overdensity, 5i, is 

6i(m,x ) = J Wpm(\x -x\)6(x, a )d 3 x. (2.17) 

Equation (2.17) is then the definition of 5i for a given (3. In our double distribution, 
(3 is free parameter, which is however constrained to be between 1 and a few on 
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physical grounds. It cannot be < 1 since the mass of the object's environment always 
includes the mass of the object itself. In fact, as (3 approaches 1, the averaging which 
produces Sg is taken only over the collapsed object itself, and inevitably returns the 
critical overdensity for collapse, <5o jC , f° r au objects. In the other extreme, j3 ^> 1, 
the average S on a scale (3m is no longer a local quantity with respect to the central 
object. When (3 grows without bound, Sg approaches for all collapsed structures, 
since averaging over the whole universe identically returns the background matter 
density, which corresponds to a vanishing density contrast. In appendix A we show 
that our double distribution becomes proportional to a Dirac delta-function around 
5e = in the limit (3 — > oo and proportional to a Dirac delta-function around 5g = 5q >c 
in the limit j3 — > 1. 

We are now ready to use random walk theory results to calculate the double 
distribution. We first find the fraction of points in space which belong to structures 
of mass between m and m + dm, which in turn are embedded in a medium of mean 
linearly extrapolated overdensity between 5g and 5g + d5g, f(m,5e,(3)d5edm. The 
double distribution then is 

_dn^_ = Pmfi f ^ ^ p)dmdS t . (2.18) 
dmdde m 

The quantity / can be written as 

fdmdS e = (f 1 dS e ){f 2 dm) (2.19) 

where fidSg is the fraction of points in space which have an average overdensity 
between Sg and Sg + d8g on a smoothing scale f3m, and f 2 dm is the fraction of points 
satisfying the previous condition which belong to collapsed structures of mass between 
m and m + dm. 

The first of the two factors above is the fraction of points still in the walk which 
are found between 5 and 5 + dS at a "time" (3m. This is the solution of the ID 
random walk problem of 5e as a function of S, with an absorbing boundary at the 
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critical collapse threshold 5o iC , as given by Eq. (2.9) but for a smoothing scale /3m, 

"(^-26 , c (q)) 2 



fidSt 



exp 



2S{f3m) 



— exp 



2S(Pm) 



d5. 



(2.20) 



The second factor (f 2 ) is the conditional probability that a point in space origi- 
nating from (f3m, 5) in the mass - overdensity plane,will reach the "wall" for the first 
time for a smoothing scale between m and m + dm. This is then the probability that 
a particular point in space is absorbed by the "wall" 5o,c(o) at a particular "time" 
S(m), provided that the origin of the walk is transferred from (0,0) to (S((3m),Si). 
This probability can then be found if, in the expression for dP(> m, a)/ 'dm (Eq. 2.11), 
we perform the substitutions £o, c — > 5o,c — 8e and S(m) — > S(m) — S((3m). Similar 
conditional probabilities were originally calculated by [2] and [3] in the context of 
rates of mergers between halos. In our case, it is 

(<5o, c (a)-<5 f ) 2 



f 2 dm 



exp 



2[S(m)-S(fim)] 



V^[S{m) - S{(3m)f /2 
Equations (2.18) and (2.19) then give 



dS 



dm 



dm . 



(2.21) 



dn 



dmdb, 



T (m,5i,f3,a) 



Pm,o 8o,c(a) —8e P 



r n ~ 




2S{/3m) 


— exp 



[S e -28 0iC (a)y 
2S(J3m) 



711 



2tt 



x exp 



[S{Pm)\V 2 [S(m) - S(f3m)} 

2 " 



3/2 



dS 



<V(a) - ^ 



2 [S(m) - S{(3m) 



dm 



(2.22) 



with 5g < 5o tC (a) and j3 > 1 so S(m) > S(j3m) 4 . Equation (2.22) is the double 
distribution we have sought and is the central result of this chapter. Integrating 
dn/(dmd5e) over 5e yields the Press-Schechter mass function, as it should. The result 
is independent of the value of (3. We explicitly perform this integration in appendix 
B. 

Note that the functional form of our double distribution is similar with that of 
the integrand used by [13] in their calculation of the cross-correlation between dark 



4 since S(m) monotonically decreases with m for all physically interesting power spectra 
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halos and mass using random walk theory, however the second variance of the field 
(corresponding to our S(f3m)) in their case refers to a fixed clustering radius and is 
independent of object mass. 

2.4 Converting Between Linear And Spherical Den- 
sity Contrasts 

In this chapter, we derive and discuss conversion relations between the linearly 
extrapolated density contrast entering the double distribution as we derived it using 
random walk theory, and the more physical density contrast predicted by the spherical 
evolution theory. In §2.4.1 and 2.4.2 we derive exact expressions for 5g(5g, a) in the case 
of the spherical evolution model, for an fl m = 1 (§2.4.1) and an Q m + Q\ = 1 (§2.4.2) 
universe (note however that all of the equations we have presented up to this point are 
cosmology-independent, and can therefore be adapted for any cosmological model). 
In §2.4.3 we present a useful approximation, valid for both types of cosmological 
models. 

2.4.1 Exact Conversion in an Einstein-de Sitter Universe 

In this section we derive a conversion relation 5o(a,5) for an Q m = 1 cosmology 
(here, 5 is the density contrast predicted for a density perturbation at cosmic epoch 
a by the spherical evolution model and 5q is the extrapolation of the density contrast 
to the present cosmic epoch using linear theory). In order to do so, we first calculate 
5(a) from the spherical evolution solution, then calculate 5q using linear theory, and 
finally require that 5(a) and 5 a (the linear-theory density contrast at epoch a) should 
agree at early times. 
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Spherical Evolution Model in an f2 m = 1 Universe 

The evolution of a spherically symmetric, overdense perturbation in an otherwise 

homogeneous Q m = 1 universe is described by the parametric equations 

2tj f sin \ 

a P = ( 12 J)2/3 ^ ~ cos6> ) ' and a = aco11 ( — 2tt — ) ' ( ^ 2 ' 23 ' ) 

where a co n is the scale factor of the universe when the perturbation formally col- 
lapses to a point, a p is the scale factor of the perturbation, and 9 is the develop- 
ment angle. Note that the perturbation will turn around (reach its maximum size, 
o p ,max = 4a co u(127r)~ 2 / 3 ) when 9 = it, at a time a = a co n/2 2 / 3 . 

The normalization of Eq. (2.23) is such that the density contrast 5 can be expressed 

as 

^fe) 3 " 1 ' (2 24) 

Hence, for any density contrast 5, Eqs. (2.23) and (2.24) can be combined to give a 
unique development angle 9(5) which is the solution to the transcendental equation 

2(1 — cos 9) 

Similarly, the spherical evolution solution for an underdensity is given by the 
parametric equations 

62/3 

a p = yl p (cosh?7 — 1) , and a = A p — — (sinhr/ — r]) 2 ' 3 . (2.26) 

where rj is the development angle in this case. Equation (2.26) together with Eq. (2.24) 

can be combined as before to give rj(5) as the solution to the transcendental equation 

6 nsinh ?7 -^/3 _ Q 
2(cosh?7 — 1) 

5o(a,5) according to the spherical evolution model 

The behavior of 5 in the linear regime in this cosmology is 

5 = 5 a. (2.28) 
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This result should coincide with the linear expansion of the spherical evolution result 
at early times. Expanding the parametric solution to second nonvanishing order in 9 
and eliminating 9, we obtain 

'12tt) 2 / 3 a 



o P (a) = a 



1 - 



20 a co n _ 

We then substite Eq. (2.29) in the definition of 5 (Eq. 2.24) to get 



(2.29) 



? 3(12vr) 2 / 3 , . 

S = \ 1 a (2.30) 
20a co ii 



which, by comparison to Eq. (2.28) gives 



^ _ 3(12tt) 2 / 3 (2 3l) 

20a co ii 

Then, the conversion relation we seek is 

fi 2 / 3 3 

h{a,5) = — [6(5) - sm9(5)] 2/3 (2.32) 

where 9(6) is given by Eq. (2.25). 

Equation (2.32) has the undesirable property that it diverges as 9 approaches 2ir. 
This is of course a consequence of the perturbation formally collapsing to a singularity 
in the spherical evolution model instead of reaching virial equilibrium. If we make 
the usual assumption that at virialization the radius of the perturbation is a max /2 
and we additionally require that 

• So(a, S) is continuous and smooth at 9 = 3tt/2 

• a p = a Pimax for all a > a coll 

then for 9 > 3n/2 (which corresponds to 5 > 9(3n + 2) 2 /8) we can replace Eq. (2.32) 
with 



8 (a, 8) = 5 ,v + 8' v (5 - 5 V ) 

3(So, c -8 , v )-(S c -8 v )(2S'^ + S' 0tC ) 



(S C - Sy) 



2 



(6 - Syf 



+ (s^y {6 - 6v) 
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(2.33) 



(see §2.4.2 for a discussion of the reasons for employing this particular functional form, 
and §2.6 for a discussion on why the effect of such a choice on the double distribution 
is negligible). In Eq. (2.33), 



6 V 



a \e=3Tv/2 
a p\e=3n/2 

a p|e»=37r/2 



- 1 



9(3tt + 2f 



18tt z 



- 1 



5' 



Ox 



S (a, <5„ 
05 



35/3 
20a 



(37T + 2) 2 / 3 



<5=<5c 



10a(l + 5 C ) 2/3 



(2.34) 



the last equality coming from the fact that after a co \\ the radius of a perturbation 
remains constant and equal to a P)max /2, while its density contrast 5 changes only due 
to the expansion of the background universe, 5 = (2a/a Pjmax ) 3 — 1 or 5 = (10a5o/3) 3 — 1. 
Finally, <5 0>c is given by Eq. (2.39) while 5' 0v is given by Eq. (2.40) for 5 = 5 V and 
9 = 3tt/2. 

To get the linear behavior of 5 for an underdensity we expand the parametric 
solution 2.26 to second nonvanishing order in rj and we eliminate rj to get 



g2/3 

a p — A P — — 



1 + 



1 a 

ToZ 



(2.35) 



Substituting Eq. (2.35) in the definition of 6 (Eq. 2.24), we get for the time dependence 
of 5 at early times, 



5 



10 A T 



(2.36) 



from which, by comparison to Eq. (2.28), we get 

3 



10A 



(2.37) 



Then, 5 (a,5) will be 



S (a,5) 



6 2 / 3 3 
20a 



[sinh?7(o") — r)(6)] 



2/3 



(2.38) 
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where rj(5) is given by Eq. (2.27). 

Equation (2.38) is valid for all 7] and its limit as 5 — > — oo is 5(5) — > —1. Thus, 
although the linearly extrapolated field can become < — 1, the corresponding value 
of the actual 5 is always > —1, as the physical requirement p p > demands. 

Critical extrapolated overdensity for collapse, So, c (a) 

The critical extrapolated overdensity for collapse can be found from Eq. (2.31) 



~ . . 3(12vr) 2 / 3 , , 
On r .{a m n) = ^ a_„ £s 1.69a 



Jo,ci«coiu = 20 " co11 ~ ^-^odi • ( 2 - 39 ) 

Note that the above equation has the functional form 5o,c( a coii) 1/ D(a co \\), where 
D(a) is the linear growth factor for this cosmology. This is also true in the f2 m +f2A = 1 
case. 

dS /dd\ a 

In addition to the relation between 5 and 5q, we will also need the derivative 
d5o/d5\ a in order to convert between true and extrapolated overdensity differentials 
in Eq. (2.13). In the case of an overdense structure, 5 > 0, Eq. (2.32) gives 



d5o 
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6 2 / 3 1- cos 9 (5) dO 

(2.40) 



10a [0(5) - sin6(5)] 1/3 dS ' 
To evaluate d9/d5 we define the auxiliary function 

F a (6, 5) = 6 2/3 (6 - sin#) 2/3 - 2(1 - cos0)(l + 5) 1/3 . (2.41) 



From Eq. (2.25) we get immediately F a (9, 5) = 0, and differentiating we get dF. t 

d6 dF a fdFA" 1 



= ^d6 + ?§fd5. Hence 



d5 85 V 09 

where 



(2.42) 



<9F a 2 1 - cos 9 

~W ~ _ 3 (1 + 5)2/3 ( 2A3 > 
23 



and 



dF» 6 2 / 3 2 1 — cos# 



- m -2(l + 5) 1/3 sm6. 



06 3 (9-sm8y/ 3 
Equation (2.40) is valid only for < 5 < 5 V . For 5 > 5 V Eq. (2.33) gives 
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5' 



0,v 



i0 3(Vc-Vv)-(^c-^)(2^ v + ^ c ) ^ 

. (*0,c + *0,v) ( J c - <*v) - 2(5 ,c - V) . . 2 

+3 {6 ~ 6v) ■ 



In the case of an underdense structure, 5 < 0, Eq. (2.38) gives 



d5o 
05 



6 2 / 3 cosh 77(0") — 1 drj 



10a [sinh 77(0") -77(0")] 1/3 <W ' 
As before, in order to evaluate dr]/d5 we define the auxiliary function 



G a (r), 5) = 6 2/3 (sinhr/ - r/) 2/3 - 2(coshr/ -!)(! + 5) 1 



/3 



Equation (2.27) implies G a (?7, 5) = so 



where 



and 



drj _ dG a /gg. 



<9G a 2 cosh 7/ — 1 
= ~3 (1 + 5) 2 / 3 

<9G a 6 2 / 3 2 coshr/-! 



(9r/ 3 (sinh 77 — rf) x / 3 



-2(l + <5) 1/3 sinhr7. 



(2.44) 



(2.45) 



(2.46) 



(2.47) 



(2.48) 



(2.49) 



(2.50) 



2.4.2 Exact Conversion in an S7 m + J7a = 1 Universe 

In this section we will derive a conversion between true and extrapolated overden- 
sity, 6(60, a) for an f2 m + f2A = 1 cosmological model. We will do so by first calculating 
the true density contrast 5(a) of a density perturbation at cosmic epoch a as predicted 
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by the spherical evolution model, then calculating 5q, which is the overdensity of the 
same spherical perturbation if extrapolated according to the linear theory until the 
present cosmic epoch, and finally requiring that at early times linear theory and the 
linear expansion of the spherical evolution model should give the same result. 

Spherical Evolution Model in an fi m + f2 A = 1 Cosmology: The Evolution 
Equation 

In the spherical evolution model, the spherical density perturbation under consid- 
eration behaves as an independent non-flat sub-universe. Its evolution is dictated by 
a Friedmann equation, 

(ift ) 2 = H2 ° Qma » (flp " 3 + W " 2) (2 - 51) 

where a p is the radius of such a spherical density perturbation in an otherwise ho- 
mogeneous universe, uj = Q\/Q m = Q" 1 — 1 (where f2 m and Q\ are the matter and 
vacuum density parameters of the background universe) and k is a constant charac- 
teristic of the amplitude and sign of the perturbation: the larger the |/c|, the larger 
the deviation from homogeneity at a given time, while a positive k corresponds to an 
overdensity and a negative k to an underdensity. Clearly then in Eq. (2.51), the first 
term in parentheses on the RHS is the matter term, the second is the vacuum term 
and the third is the curvature term, which can have a positive or negative sign depend- 
ing on whether we are studying an "open" (underdensity) or "closed" (overdensity) 
perturbation. The normalization of a p is such that, had the specific spherical region 
begun its evolution with no curvature (k = 0), a p at the present cosmic epoch would 
have been a p ( K= o),o = 1- For this reason, the density contrast 6 of the perturbation at 
epoch a is given by Eq. (2.24) 

The behavior of the perturbation radius Clp cLS cl function of the universe scale factor 
a can be found by taking the ratio of the Friedmann equations of the perturbation 
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and the background universe, thus obtaining [66] 



da p 
da 



a n 1 + ioa^ — k a wai — fta D + 1 



ua z 



(2.52) 



Equation (2.52) implies that the smallest positive perturbation which will turn around 
and collapse corresponds to the smallest positive k for which the equation 



CJOp — Kdp + 1=0 



has a real positive solution [9] . This gives 



^min,coll 



3^ 1 / 3 /2 2 / 3 . 



(2.53) 



(2.54) 



Equation (2.52) can then be re-written as 



da p 
da 



-l , 2 \ 1/2 



K < K m in,coll or 

K ^ K min,colh a < Ota 



1 9 



1/2 



K ^ ^min.coll; a > °ta 



(2.55) 



where at a is the scale factor of the universe when the perturbation reaches its maxi- 
mum (or turnaround) radius. The turnaround radius is the smallest of the two positive 
solutions of Eq. (2.53), 



7T 



(2.56) 



Equation (2.56) has an asymptotic behavior a p>ta ~ l/« when k/uj ' ^> 1, as expected 
from Eq. (2.53). The maximum possible turnaround radius, a Pj ta,max is achieved for 



k — ^min coil and is a 



p,ta,max 



(2a;) 1 ' 3 . All other collapsing overdensities will have 



^p,ta ^ *^p,ta,max- 
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Qualitative Description of the Evolution of Structures 

The introduction of the additional vacuum term in the Friedmann equation consid- 
erably complicates the simple classification of density perturbations to overdensities 
(all of which turn around and collapse in an Q m = 1 cosmology) and underdensities 
(all of which expand forever). In the Q m + = 1 universe there exist overdensities 
which will continue to expand forever. The behavior of a perturbation in such a 
cosmology is parametrized by the quantity k/lu 1 ^ 3 , and we can identify the following 
cases. 

Case I, k/uj 1 / 3 < —1: large underdensities, expanding forever. The table 
below shows the relative magnitude of the three terms in the Friedmann equation 
(matter, curvature and vacuum term) for different values of the scale factor of the 
perturbation. The first line in the table indicates the hierarchy of the three terms, 
from largest to smaller, for each range of the scale factor. The second line indicates 
the dominant term in each scale factor range. The third line shows the approximate 
dependence of the radius of the perturbation, a p , on time, assuming that only the 
dominant term contributes to the Friedmann equation in each range. 



°P < \k\ 


|«| <a p < w 


i < a < /H 


^ fVA 


MCV 


CMV 


CVM 


VCM 


matter 


curvature 


curvature 


vacuum 




a P ~ t 


a p ~ t 


a p ~ e* 



Case II, —1 < k/cj 1//3 < 1: small perturbations, expanding forever. 

These can be either underdensities (k < 0) or overdensities (k > 0). In both cases 
the curvature term never becomes dominant. The following table shows their different 
evolutionary stages (as in Case I). 
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a D < \ P-*- 
f y w 


/N < a < i 


< a D < ji, 

V w p M 


a P > |k| 


MCV 


MVC 


VMC 


VCM 


matter 


matter 


vacuum 


vacuum 


a p ~ t 2 / 3 


a p ~ t 2 / 3 


a p ~ e* 


a p ~ e* 



Case III, 1 < k/cj 1 / 3 < 3/2 2//3 : "coasting" overdensities, expanding 
forever. These overdensities continue to expand forever despite the fact that they 
go through a phase in their evolution when the curvature term becomes dominant 
and their expansion slows down. During this phase, the contributions of the matter 
and vacuum terms, which are the ones driving the expansion, add up to a value 
always higher than the curvature term, although the curvature term is larger than 
each one of them. When the perturbation enters the curvature-dominated phase, the 
expansion rate decreases and the perturbation grows much more mildly than t 2 ^ 3 . 
The expansion rate reaches a minimum at a p = (2cj) -1 / 3 , after which it increases 
again as the perturbation approaches the phase of exponential expansion. This phase 
between the matter- like expansion and the exponential expansion is denoted by (*) 
in the table below. 



a D < - 

V K 


i < a p < 4= 

ft f yuj 


-S7= < On < \/^ 

yui v V u 


a D > \/* 

v v u> 


MCV 


CMV 


CVM 


VCM 


matter 


curvature 


curvature 


vacuum 


a p ~ t 2 / 3 


(*) 


(*) 


a p ~ e* 



Cases I-III are all sub-cases of the Lemaitre model ([67], [68]), which features an 
inflection point at a PiC = (2u;)~ 1 / 3 where a p = while d p > 0. The rate of expansion 
initially decreases to achieve its minimum (positive) value when a p = a P)C , after which 
point the expansion accelerates again. 

Special Case, k/oj 1//3 = 3/2 2,/3 : Eddington Overdensity. This overdensity 
is the lowest k overdensity which does not expand to an infinite radius. However, it 
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does not turn around and collapse, but it approaches its (finite) turnaround radius, 
<2p,max = (2cj) -1 / 3 (from Eq. 2.56) as t — > oo. As seen by an observer inside this 
overdensity, as t — > oo the part of the universe outside a Pjmax will accelerate away and 
eventually exit the horizon, and the observable universe ("local Eddington bubble") 
will asymptotically approach the Einstein static universe (as in the Eddington model 
with a cosmological constant). 

Case IV, k/uj 1 / 3 > 3/2 2//3 : large overdensities, eventually collapsing. 
When such a structure enters the dominant-curvature-term phase, its expansion rate 
starts to decrease (a p ~ t e with e = e(t) monotonically decreasing from 2/3 to 0), until 
the expansion halts, at a p = a pta which occurs at a time t ta , given in table 2.1. After 
t ta the perturbation turns around and contracts, its evolution being symmetrical in 
time about t ta , i.e. a p (t) = a p (2t ta — t) for t > £ ta (this is a consequence of Eq. (2.51) 
and holds for any cosmological model as long as the RHS of the Friedmann equation 
involves no explicit time-dependence). Eventually, the perturbation will formally 
collapse to a singularity at time t co \\ = 2t ta - 



«p < i 


- < a p < a Pi ta 






MCV 


CMV 


CMV 


MCV 


matter 


curvature 


curvature 


matter 


a p ~ t 2 / 3 


a p ~ t e 


a P ~ (2t ta - t) e 


a p ~(2t ta - t) 2 / 3 


expansion 


expansion 


contraction 


contraction 



In all of the cases discussed above, the transitions between different phases of 
their evolution occur at characteristic times, those of matter- vacuum equality £mvj 
matter-curvature equality £mc and curvature- vacuum equality tcv- At these times 
(shown in table 2.1), the corresponding terms in the Friedmann equation become 
equally important. Note that in the case of the Eddington overdensity and of case IV 
collapsing overdensities, matter-vacuum equality and curvature-vacuum equality are 
never reached, and the vacuum term never dominates over any of the other terms. 
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Table 2.1: Characteristic times of the spherical evolution model in an fi m + = 1 
universe. 

event time 

turnaround t ta = „ \„ f n ap ' ta da 



H \/n^ JO PV Loa^-Kap+l 



matter- vacuum 
equality 



1 r^" 1/3 a 

''MV - w „./o- Jo "°F 



H Vn^ JO <~"vy 



matter-curvature 
equality 



* MC - HoVifc So ' rfa P J ^-laTTT 



curvature-vacuum i , 



equality 



In the next section we derive exact solutions for the time-evolution of a p for 
perturbations of different curvature. However, surprisingly accurate approximate 
solutions can be derived using only linear theory and Eq. (2.81). Solving for the 
spherical collapse density contrast we get 

.r 

Since a p = a(l + 5 a ) -1 ' 3 , we can write for collapsing overdensities 

Sc/3 



S a « I 1 - ^ I - 1 . (2.57) 



a p ~ a 



1 (5 c D(a)/D(a c 



5 C /3 



(2.58) 



^(Ocou). 

where the initial conditions (curvature) of the perturbation are parametrized by its 
collapse epoch, et co u, while the cosmology enters through the functional form of the lin- 
ear growth factor and the linear collapse overdensity, 5 C . Similarly, for perturbations 
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which expand forever we can write 



a p ~ a 



S D(a)/D(a ) 



Sn D(a) 



6c/ 3 



-i <5 c /3 



(2.59) 



5 C D(a ) 

where the curvature of the perturbation is parametrized by its extrapolated linear 
density contrast at the present epoch, 5q. Note that for overdensities which expand 
forever, 5q > and a p < a, while for underdensities 5q < and a p > a. Also, because 
D(a) asymptotes to a constant value for a — > oo (as we will see in the next sections), 
a p grows proportionally to a at late times. This is the exponential expansion phase, 
described in our analysis above. 



Solutions of the Evolution Equation 

For eventually collapsing structures (ac > « m m,coii) > separation of variables in 
Eq. (2.55) and integration yields, 



o J uy 3 + 1 



\/u)X :i — KX + l 



a < Ota 



(2.60) 



o r ^Vxdx r °P Vx-dx a>a 

JqJujx 3 — kx+1 Jo\/wx 6 — kx+1 — 



where a ta is the cosmic epoch when a p = a Pj t a - Now the integral on the LHS of 
Eq. (2.60) can be calculated using [9] 



Vydy 



Jury* + 1 3 
The integral of the RHS can be re-written as 



-u x l 2 sinh 1 \/i 



uy= 



(2.61) 



"p 



yfxdx 



(2.62) 



'o V UJX 3 — KX + 1 3 

where Vi is the incomplete vacuum integral of the first kind, defined in appendix C, 
and 

r = a p /a pM , fj, = (waJJ -1 • (2-63) 
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Note that for this range of curvature values, k/uj 1 ^ 3 (which is the quantity which 
parametrizes the behavior of the perturbation with time) is a function of fi alone, 
with 

k/u 1 ' 3 = (l + /i)//i 2 / 3 . (2.64) 
Using Eqs. (2.61) and (2.62), Eq. (2.60) can be rewritten as, 

^^{^[^(r,/!)]} 2 / 3 , a<a ta 

(2.65) 

^{sinhpViCl.^-ViCr./i)]} 2 / 3 , a>a ta 

Equation (2.65) is the spherical evolution solution for a collapsing perturbation with 
amplitude k (parametrized above by /i) in an ^a ~t~ — 1 •> ^A/^m — ^ universe, 
and gives r (and hence a p ) as a function of a for this model. Note that Vi(r, //) is the 
development angle for this cosmology. 

The scale factor of the universe at turnaround for a given collapsing overdensity 
can be found immediately from Eq. (2.65), 

a ta = c^ 1 / 3 [sinh Vx (1, /i)f 73 . (2.66) 

The scale factor of the universe at collapse, a co n (when the scale factor of the perturba- 
tion becomes formally zero, a PjC = r c = 0) is, from Eq. (2.65) and since Vi(0,/i) = 0, 

Gcoii = ^ 1/3 [sinh2V!(l,/i)] 2/3 . (2.67) 

Equation (2.65) should not be applied "literally" until the final collapse of the 
perturbation to a singularity, since the physical picture for the late stages of the 
evolution of a perturbation involves virialization at a finite radius. It has been shown 
by [69] that the analogous arguments which give a PiV = a Ptta /2 for the Q m = 1 
cosmology give, for an Q m + = 1 universe, 

, 2 + 2uoal ta 
Aual v —a v ^ + 1 = 0. (2.68) 



p,ta 
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The physically meaningful solution of 2.68 which gives the correct behavior for uj — > 
is (using Eq. 2.63) 



2/^ + 2 1 _j / 27/^ 



(2.69) 



The scale factor a v of the universe when the scale factor of the perturbation past its 
turnaroundbecom.es equal to a p>v , will be given by the second branch of Eq. (2.65), for 
a p = a P ,v Then, the validity range for the second branch of Eq. (2.65) is a ta < a < a v . 

For a > a v , we can no longer use the spherical evolution solution to describe the 
physical picture of interest (virialization). In the next section we will present a simple 
recipe we will use to follow the late stages of evolution of the perturbation which 
satisfies the desired boundary conditions (a p = a pv at a co \\ and constant thereafter). 

For perpetually expanding structures (k < K m i n , C oii)> Eq. (2.55) gives 

na Vy d y _ r V^dx 



\/ooy 3 + 1 Jo \fujx 3 — kx + 1 
In this case, the integral on the RHS can be rewritten as 

r , f d \, = V i/2 ^, ^ (2.7i) 

Jo y/ojx* — kx + 1 o 
where TCi is the hyperbolic vacuum integral of the first kind, defined in appendix C, 
and 

r = a p /\a PtR \ , w = (^Ig^rJ 3 )- 1 , (2.72) 

where a p R is the only real (and always negative) root of Eq. (2.53) when k < /t m i n;CO ii, 

i 

—LJ3 

g p ,R = i = j = (2-73) 



As in the case of collapsing perturbations, k/uj 1 ^ 3 is a function of w alone, with 

k/uj 1/3 = (1 - w)/w 2/3 . (2.74) 
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Then, w = 1 is a flat subuniverse (not perturbed with respect to the background), 
and perturbations with w > 1 are underdensities while 1/4 < w < 1 correspond to 
non-collapsing overdensities. 
Then, Eq. (2.70) becomes 

a = u~ 1/3 {smh[Hi(r,w)]} 2/3 (2.75) 

which is the spherical evolution solution for a non-collapsing perturbation with am- 
plitude a (parametrized above by tu) in an Q\ ~\~ fi m — 1, fi^Z^m 

= uj universe. As 

for collapsing perturbations, 7i\ (r, w) is the development angle. 

<5o(<5) according to the spherical evolution model 

The linear theory result for a growing-mode perturbation in an Q m + Q\ = 1 
cosmology is [70] 

where D, the linear growth factor, is given by D(a) = A[(2a;) 1 / 3 a] with 

To find the relation between 5 an d we expand the exact (5 = a 3 /a 3 — 1) and linear 
relations for the overdensity to first order in a and demand that the coefficients be 
equal. Thus we get [9] 

\ d = o A M/sA ' I 2 ' 78 ) 



3A[(2cu) 1 /3 ao ] u 3A[(2uy/ 3 } 

since ao = 1 . This is the value of the constant k for a perturbation which at a 
cosmic epoch ao = 1 has a linearly extrapolated overdensity 5q- Note that since the 
linear theory result is the same for both underdensities and overdensities, Eq. (2.78) 
holds for both cases. For underdensities, both k and 5 will be negative, while for 
overdensities both will be positive. 
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Table 2.2: Applicability limits for different branches of the spherical evolution solu- 
tion, where: c5 Ed (a) is the density contrast of the Eddington overdensity; c5 ta (a) is the 
density contrast of an overdensity turning around at a; 5 v (a) is the density contrast of 
overdensity reaching its virial size at a; and S c (a) is the density contrast of virialized 
overdensity formally collapsing to a point at a. 



Limit 


Expression 




Auxiliary Relations 




5 E d(a) 


fed (a) = 


2uj( a , y i 

\ r Ed(a) J 


sinh -1 y/uoa 3 — Vi(r E d, 2) 


= 


<5ta(a) 


<5ta( a ) = 


ua 3 fi tau (a] 


~ 1 


sinh" 1 y/ua 3 - Vi(l,/i ta ) 


= 


6y(a) 


6 v (a) = 


( a 


) 3 -l 


sinh -1 \/cua 3 — 2Vi(l,/x v ) 
r v (/i) = a PiV (/i)/a Pi t a (/i) 


+ Vi [r v (/i v ),/i v ] = 


\ a p,v[Mv(a)] 


8 c (a) 


S c (a) = 


a 


3 

- 1 


sinh -1 y uja 3 — 2Vi(l,/x c ) 


= 





At any given epoch a, there is a unique perturbation (parametrized by k or, 
equivalently, by \x or vo) which will have achieved a true density contrast 8 at that 
time. Therefore, to calculate the desired conversion relation 5o(a, S) we first calculate 
k (or, equivalently, \x or w) from the given a and 5 and the appropriate solution of 
the evolution equation (2.65 or 2.75). Then, we use Eq. (2.78) to evaluate Sq. 

To determine which is the appropriate solution of the evolution equation we need 
to use for each 5, we calculate the limits of applicability of each equation in terms of 
S. 

• Equation (2.75) is applicable for all forever expanding perturbations. For any 
given epoch a, the maximum density contrast of such perturbations is achieved 
by the Eddington perturbation and is equal to <5 E d( a ) ; given in line 1 of Table 
2.2. Then, the applicability domain of Eq. (2.75) is — 1 < 5 < 5^(0), and the 
conversion relation in this case takes the form shown in the 1st line of Table 
2.3. 

• The first branch of Eq. (2.65) is applicable for eventually collapsing pertur- 
bations which, however, have not reached their turnaround radius yet. The 
maximum 5 of all such perturbations at a given a is achieved by the pertur- 
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bation which is turning around at a, and is equal to S ts .(a), given in line 2 of 
Table 2.2. The applicability domain of the first branch of Eq. (2.65) is then 
^Ed(a) < S < 5ta(fl) and the conversion relation in this case is shown in the 2nd 
line of table 2.3. 

The second branch of Eq. (2.65) is applicable for eventually collapsing pertur- 
bations which are past their turnaround but which have not yet reached their 
virial radius. The maximum 6 of such perturbations at a given a is achieved 
by the perturbation which is reaching its virial size at a, and is equal to S v (a), 
given in line 3 of Table 2.2. The applicability domain of the second branch of 
Eq. (2.65) is then <5ta(o) < 5 < 5 v (a) and the conversion relation in this case is 
shown in the 3rd line of table 2.3. 

Perturbations which have reached their virial size but have not yet reached 
their designated collapse time, a co n(/i), need to be treated separately, since the 
spherical collapse model fails (does not agree with the physical picture we would 
like to describe, although it is still formally applicable) for radii smaller than 
the virial radius. Since a realistic description of the microphysical dissipation 
processes which lead to virialization is far beyond the scope of this analyti- 
cal calculation, we will adopt a prescription which is driven by mathematical 
simplicity. We will assume that for 5 v (a) < 5 < 5 c (a) the conversion relation 
5o(a, 5) has the simplest polynomial form which satisfies the following physically 
motivated boundary conditions: 

— The extrapolated overdensity is continuous and smooth at 5 V , so 5o{a, S v ) = 
5o tV and 85^/35 = 5' 0v as given by the appropriate relations of the pre- 
vious branch (3rd line of tables 2.3 and 2.4 correspondingly). 

— After time a co n(/i) the radius of the perturbation remains constant and 
equal to the virial radius so changes in the (true) overdensity are only 
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due to the increase of the scale factor of the background universe. This 
then implies that 5o(a,5 c ) = 5o, c given by Eq. (2.79) and d5o/dS = 
(d5 /dfi\ llc )(dfi/d5\s c ) = 5' 0c given in Table 2.3 line 4. 

The conversion relation in this case is shown in Table 2.3 line 4, while its appli- 
cability domain is 5 v (a) < 6 < S c (a), with 5 c (a) given in Table 2.2 line 4. Past 
their collapse time, perturbations are treated as virialized objects without sub- 
structure and are not relevant as "local environment" of other objects for the 
purposes of our double distribution calculation, hence it is not necessary to have 
a conversion relation of 5 > 5 c (a). The calculation would be simplified (this last 
branch would be unnecessary) if we chose to regard perturbations as virialized 
objects after the moment they reached their virial size after turnaround, at time 
a v . However, we will retain the usual assumption that objects virialize at time 
a co n for consistency with existing Press-Schechter calculations. 



Critical extrapolated overdensity for collapse, S , c (a) 

We need to find the critical So iC (a) for collapse if the field is linearly extrapolated to 
the present epoch, i.e. the value the linearly extrapolated to the present overdensity 
must have, for a structure to have collapsed at universe scale factor a. This, from 
Eq. (2.78), will be 

6 (a)- 3A ^) 1/3 ] 1 +^H (2J9) 

where again fi c (a) is given by Table 2.2 line 4. 

The dependence of 5o, c on a can also be expressed in terms of the linear growth 
factor, D(a), as was the case for the Q m = 1 universe. The conversion relation between 
5(a) and 5 a (the linear-theory result for the density contrast at time a) is independent 
of a. In other words, as long as 6 and S both refer to the same time, knowledge of the 
one uniquely defines the other, independently of the actual time at which they are 
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Table 2.3: Different branches of conversion relation 5o(a,5) for an Q m + Q\ = 1 
universe, where: branch I is — 1 < S < 5va\ branch II is 5™ < S < <5 ta ; branch III is 



>Ed 

<5ta < 5 < <^vj an d branch IV is S v < S < S c . 



Branch 5q(5, a) 



Auxiliary relations 



I 


3-A[(2w) 1/3 


1— za(a,S) 


sinh 1 V tua 3 


-n x 


a 




= 


2V3 


[ro(a,5)l 2/3 


II 


3A[(2w) 1/S 


1+Ai(a,5) 


sinh -1 y/ua 3 


-Vi 


a( 




= 


2V3 


MaM 2/3 



III 



3A[(2aOV3] i+^( a ,g) 
2 1 / 3 [ M (a,5)] 2 / 3 



sinh -1 Vwa 3 " - 2Vi(l,/i) 



+Vi 



(^) 1/3 ^ 







IV 



<^o,v = ^o( a ) <^v) (this Table line 3) 
(Table 2.4 line 3) 



7/ _ 95 
°0,v - g5 



(5c-5 v ) 

x(5-5 v ) 2 

(^,c+^o.v)(^-M-2(5 ,c-<5o,v) 



V = 5 (a, <f c ) (Eq. 2.79) 



0,c 



95o 
95 



(<5c-<5v) 3 

x(5-5 v ) 3 



5 C 



(2w)V3] {6^al^ 2 J s -l-^ 



(2 M 2) 1 / 3 (5 C +1) 

/x c (a) in Table 2.2 line 4 



both evaluated. As 5(a) — > 00, 5 a — > 5 C , the linear-theory density contrast at the time 
of collapse (given by Eq. 2.79 for a = 1). Therefore S c is the same for perturbations 
of all curvatures, and, using Eq. (2.76), we can write 

5 >coii) = ■ (2.80) 

dS /dS\ a 

In addition to the relation between S and So, we will also need the derivative 
dSo/dS\ a in order to convert between true and extrapolated overdensity differentials 
in Eq. (2.13). The calculation is similar as in the case of an Q m = 1 universe, and the 
results are summarized in table 2.4. 

In table 2.4, Hiir^w) is the hyperbolic vacuum integral of the second kind, and 
V2(r, n) is the incomplete vacuum integral of the second kind, defined in appendix C. 
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Table 2.4: Different branches of derivative dSo/dS\ a for an Q m + Q\ = 1 universe. 
The roman numerals correspond to different branches as in Table 2.3 



Aux. Function ^ 



(a, 6) = Additional Functions 



$x = sinh 1 \fuocfi 
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sinh 1 V ua 3 —Hi 
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$ 2 = sinh 1 \fuja? 
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95 ^flfo/ 
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,3/2 



2(1+5)3/2 



a 3 a(M+l) 

(1 + «T M 2/3(l + «)l/3 



(2 M 2 ) 1 / 3 



'3A[(2w)V3] K 2 
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sinh" 1 vW-Vi [a {-^Y^fi 



$ 3 = sinh 1 \fua? 
HI -2Vi(l,|i) 
+V(r,/x) 
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o / a 3 "(m + 1) I 1 

-*'(!+«) ^2/3(1+4)1/3 ^u, 



9*3 = (2m 2 ) 1/3 
95„ 3A[(2o;) 1 /3] 



v 2 



«(S)Va» 



6m rfVi(l.M) ] 
/i— 2 d/i J 

sinh -1 vW-2Vi(l,/i) 



+Vi 



?/ | O 3(^o,c-5o,v)-(5c-<5v)(25 ( ') v +<5o, c ) / r r\ 
V + ^ l d - d vJ 

i o ( < 5o.c+^o,v)( <5 c-5v)-2(5() jC -(5() iV ) , , 2 

jy +d (ST^M 3 1<>-<M 

with 5 ,c V <% ;C , ^ 
as in Table 2.3 line 4 
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2.4.3 An Approximate Conversion Relation 

An excellent approximation to these conversion relations can be derived from the 
expression 

5 a ^5 c [l-a + .y- 1 ^] . (2.81) 

Similar approximations were suggested by [71] and [15]. Equation (2.81) relates the 
linear overdensity at a time a to the true overdensity at the same time, and its 
accuracy is better than 2% throughout its domain for both Q m = 1 and Q m + Q\ = 1 
cosmologies. Its functional form is much simpler and more intuitive than the more 
accurate fit of [13]. The cosmological model enters only through 5 C . For the Einstein- 
deSitter universe, 5 C is given by Eq. (2.39) for a co n = 1, while for the Q m + = 1 
universe it is given by Eq. (2.79) for a — 1. Note that 5 C is related to the quantity 
<W( a ) (which appears explicitly in our double distribution expression) through 

*o» = (2.82) 
D(a) 

where D(a) is the linear growth factor in the relevant cosmology. 

The limits of Eq. (2.81) are the same as the ones required for the exact conversion 

relation. When \5 a \ <C 1, 6 a ~ 5 a . In addition, S a — > — oo as S a — > — 1, and 5 a — > oo 

as 5 a — >■ Re- 
using Eq. (2.81), 

~5 e ^^j4~S c \l-(l + S £ )- 1/l ] , (2-83) 
D{a) L J 

where a is the time at which we want to evaluate the double distribution. Note that 
close to virialization, Eq. (2.83) loses its applicability (as does the spherical collapse 
model), and has to be replaced by a recipe which does not diverge in 5^. We have 
presented such recipes in §2.4.1 and 2.4.2, however the exact functional form of the 
conversion relation in this regime cannot affect any of the physically interesting results 
as the amplitude of the double distribution decreases rapidly enough with 5g that the 
contribution of the high-delta tail to the integrated mass function is negligible. We 
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have verified this fact by comparing the integral of our double distribution over Se 
with the Press-Schechter mass function. When we extended the integration up to <$o )C , 
the results agreed to the accuracy of the numerical integration. When we extended 
our integration only up to 5o,v (just below the application of our virialization recipe), 
the error relative to the Press-Schechter mass function was less than 0.02%. 

2.4.4 Clustering Scale Lengths and Correction for Central 
Object Contamination 

The definition of (3 and 5e described above was sufficient for us to derive the 
double distribution from random walk theory. However, from a physical point of 
view, the presence of a collapsed structure at the center of the "environment sphere" 
contaminates the evaluation of the average "environmental" overdensity. If we want 
the double distribution to describe the properties of matter surrounding collapsed 
objects, we need to correct for the presence of the objects themselves. 

We will employ a simple, "top-hat" physical picture to calculate an appropriate 
correction (see also [15]). Note however that our correction is approximate, since the 
filter we used to smooth the overdensity field was k— sharp rather than top-hat in 
space. 

Let S c be the (true) overdensity of a collapsed object of mass m and radius _R V , and 
Se be the overdensity of the "environment sphere" of radius R c . The "environment 
sphere" encompasses a mass /3m, including the central collapsed object. We want to 
find the average overdensity 5 cxt of that part of the "environment sphere" which is 
external to the central object. For the collapsed object we can write 

m=^7iR 3 v (l + 5 c )p m , (2.84) 

where p m is the mean matter density of the universe at the epoch of interest. For the 
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environment sphere, including the central object, we can write 

Pm = ^7iR 3 c (l+S e )p m . (2.85) 

From Eqs. (2.84) and (2.85) we get fi* = R 3 C {1 + 5 e )//3{l + 5 C ). It thus follows that 
the length scale R e associated with the clustering parameter (3 is 



i + s t 

3(3m 



(2.86) 



,47r(l + ^)p n 

We see that for a fixed clustering parameter (3, the length scale associated with 
an object of mass m is mass-dependent, scaling linearly with the virial radius but 

1 /3 

larger by a factor [f3{l + <J C )/(1 + 5 e )] 1/0 > 1. Thus, we can roughly think of the 
clustering scale parameter as a measure of how many virial radii we include as the 
local environment around each structure 5 . 

Having identified the environmental length scale, we can now isolate the environ- 
mental overdensity from that of the collapsed object. The volume of the environment 
sphere external to the central object contains a mass 

{P - l)m = ~7r(i^ - Rl){\ + S cxt ) Pm . (2.87) 

Using Eq. (2.86) to eliminate R v , and dividing by Eq. (2.85) we obtain 

5cxt " 0(l + 6 c )-(l + 6 t ) - 1 ' (2 ' 88) 
which is the contamination-corrected overdensity for an environment sphere with 
uncorrected overdensity 5g. Then, the contamination-corrected double distribution 
will be given by 

dn . dn dSg 

-(dcxt, m, a) = - — — [^(dext, a),m, a\ — — , (2.89) 



amadext amade do cxt 



5 Note however that for fixed (3, the number of virial radii included in the environment depends 
on Si and is larger for underdense environments. 
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where 



x/x ^ /3(l + <U)(l + 4) , , 9Qm 



and 



dS cxt [(p -1)(1 + 5 C ) + + 



2 ' 



(2.91) 



2.5 Derivative Quantities of the Double Distribu- 
tion 

We now have enough tools to calculate derivative quantities of interest. The 
number density of collapsed objects of mass greater than some minimum m m i n 6 
embedded in a medium of local overdensity between <5 ext and 5 cxt + d5 cxt is 

"T7 (> "lmin)aclcxt = CtClcxt^^T / ~ "-"Om , (2.92) 

dOext ooe ao cxt J m=mrnii fLmdo£ 

while the density of matter in collapsed objects of mass > m min embedded in a 
medium of local overdensity between <5 ext and 5 0xt + d5 cxt is 

dp 88 1 d5i f°° dn , 

— — (> Wmm)QOext — d8 ext / m- — -^dm . (2.93) 

"0 ex t C0£ CXOext Jm=m min dmaO£ 

Of all the matter in the universe which belongs to collapsed objects of mass > 
m min , the fraction by mass which lives in underdense neighborhoods is 

f _ i"L=-l m min)^cxt 
J 5ext= _l d 5cxt l> ^minjaCext 



6 The introduction of a finite minimum mass m m i n is necessary for both physical and technical 
reasons. Physically, the mass of collapsed objects is strictly forced to have a lower bound, not only 
due to the finite mass of the dark matter particle, but also due to the existence of a dark matter 
Jeans mass, however small this may be. In addition, the dark matter perturbation transfer function 
imposes cutoffs at scales of order of an earth mass [????]. Practically, the Press-Schechter 
dn/dm diverges as m — ► and setting a minimum mass is required to extract interesting information. 
For the purposes of this thesis, the selected mass cutoff will generally be such that the population 
of interest will include super-galactic scales only. 
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Then, the mass fraction of the matter defined above which lives in overdensities will 

be fp,ov 1 fp,un- 

Similarly, of all the objects with mass m > m m ; n , a fraction by number which lives 
inside underdensities is 

, _ /alt=-i m min)^cxt 

L=-id£(> m -)4xt 

The complementary number fraction of such structures living inside overdensities will 

be fn,ov 1 fn,un- 

The number-density-weighted mean 5 cxt for structures of mass > m min is 



IsL=-l dfc^ m min)^ext 



with a variance 



2 sL=-i^ - ^)^) 2 dfc(> "w)^«t 

5 ' n = ~i \7x ' ( ' 

Similarly, the matter-density-weighted mean 5 for structures of mass > m min is 

W- 9 = — ~t ^ ^ > 



with a variance 



//l=-idfc( >m --)^-t 

2.6 Results 

In this section we present plots of the double distribution itself as well as of 
its various physically interesting derivative quantities. We compare results derived 
for a concordance, Q m + Qa = 1 universe with WMAP parameters (as = 0.84, 
h = 0.71, tt m = 0.27, tt h = 0.04, [24]), and for an Einstein-deSitter (Q m = 1) 
universe with h = 0.71, fit = 0.04, but a$ = 0.45. The different power-spectrum 
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Figure 2.1: Surface plots of the double distribution for z = and /? = 2 in f2 m +^A = 
1 (left panel) and Einstein-deSitter (right panel) universes. The mass is measured in 
M . The vertical axis is linear, with the m— 5 axes level corresponding to dn/dmd8i = 
and the highest point corresponding to dn/dmd5e = 2.72 x 1CT 2 (left panel) and 
dn/dmd5i = 1.89 x 10" 1 (right panel) objects per Mpc 3 per 1O 15 M . 

normalization in the Einstein-deSitter case was selected so that the Press-Schechter 
mass function in this case coincides with that of the concordance universe on a mass 
scale of 5.5 x 1O 14 M , which is between the values of m 8 (mass included in a sphere 
of comoving radius 8h~ x Mpc) for the two cosmologies (mg = 2 x 1O 14 M for the 
concordance universe while mg = 8 x 1O 14 M for the Einstein-de Sitter universe). 
This value of <jg is also consistent with the fits of [9] given the WMAP result for the 
concordance universe. Finally, we use fitting formulae of [72] for the adiabatic cold 
dark matter transfer function to calculate the density field variance S(m). In this 
section, 5 always refers to <5 ext , the true overdensity of that part of the "environment 
sphere" which is external to the central object. 

Figure 2.1 shows a 3-dimensional rendering of our double distribution as a function 
of mass and overdensity for fixed (3 = 2 and z = 0. The left panel corresponds to 
the concordance universe while the right panel corresponds to the Einstein-deSitter 
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universe, and this arrangement is retained throughout this section. 

The features of the double distribution are demonstrated in more quantitative 
detail in Figs. 2.2-2.5. Figures 2.2 and 2.3 show slices of the double distribution at 
constant values of mass. In Fig. 2.2, different curves correspond to different values 
of the central object mass. In Fig. 2.3, all curves are for an object mass m = 5.5 x 
1O 15 M , and different curves correspond to different redshifts. Their most prominent 
feature is the pronounced peak at a relatively low value of \8\, indicating that for 
each given pair of z and m, there is a preferred, "most probable" value of the local 
environment density contrast. As we can see in Fig. 2. 2, the location of this peak 
moves to higher values of the density contrast as the mass of the object increases: 
small structures are preferentially located in relative isolation, while larger structures 
are more likely to be found in clustered environments. This result fits well in the 
picture of hierarchical structure formation, as smaller structures tend to be merged 
into higher-mass objects as time progresses. Lower-mass objects which are initially 
part of underdensities are less probable to undergo mergers, and hence are more 
likely to survive at late times than objects which are initially part of overdensities. 
Conversely, higher-mass structures are more likely to be parts of overdensities where 
they can accumulate mass more easily through mergers with smaller structures. 

Note, however, that in the hierarchical structure formation picture, the mass scale 
where the exponential suppression of collapsed structures sets in increases with time. 
Thus, any given mass scale starts out as being a "high mass" at early times and even- 
tually becomes a "lower mass" as it enters the power-law regime of the Press-Schechter 
mass function. Hence, according to the argument we used to explain Fig. 2.2, the 
double distribution for any given mass scale should peak at increasing S values with 
increasing redshift. This is because a particular mass scale used to be closer to the 
high-mass end of the halo distribution in the past than it is today. Indeed, this is 
the trend seen in Fig. 2.3. As we would expect, the peak of the distribution moves 
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Figure 2.2: Slices of the double distribution function at various fixed values of the 
mass for f2 m + f2 A = 1 (left panel) and Einstein-deSitter (right panel) universes. The 
units of the double distribution are number of objects per Mpc 3 per 10 15 M Q . 




Figure 2.3: Slices of the double distribution function at m — 5.5 x 10 14 , M and 
for various values of redshift z, for Q m + £1\ — 1 (left panel) and Einstein-deSitter 
(right panel) universes. The units of the double distribution are number of objects 
per Mpc 3 per 1O 15 M . 
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Figure 2.4: Slices of the double distribution function at constant values of 5 for 
z = 0, (3 = 2 and for Q m + = 1 (left panel) and Einstein-deSitter (right panel) 
universes. Solid line: 5 = —0.5; dashed line: 5 = 0; dot-dashed line: 5 = 0.5; double- 
dot-dashed line: 5 = 3. The units of the double distribution are number of objects 
per Mpc 3 per 1O 15 M . 

to higher 5 values with increasing redshift. The significantly more pronounced sup- 
pression of this mass scale in high redshifts in the Einstein-deSitter universe is due to 
the different power-spectrum normalization in the two cosmological models. Because 
of our choice in the power-spectrum normalization, the exponential suppression in 
the number density of structures sets in at low masses in the Einstein-de Sitter case 
than in the concordance universe. Thus, there is a tendency to see more structures 
of higher mass in our concordance results than in the Einstein-de Sitter case, de- 
spite the intuitive expectation that a higher Q m universe should have more massive 
structures at late times due to its ability to continue to form structures even at the 
present epoch. This would indeed have been the case if the power-spectrum had been 
normalized in the same way. 

That halos of a given mass are more strongly clustered with increasing redshift 
was also found by [73], who used Ag(m) (the rms overdensity in the number of haloes 
more massive than some mass scale after smoothing with a spherical top-hat filter of 
comoving radius 8/i _1 Mpc) as a measure for halo clustering. A tendency of higher 
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Figure 2.5: Slices of the double distribution function at constant values of S for 
z — 0, (3 — 2, and for Q m + Q\ = 1 (left panel) and Einstein-deSitter (right panel) 
universes. Solid line: 5 = 10; dashed line: 5 = 20; dot-dashed line: 5 = 30. The units 
of the double distribution are number of objects per Mpc 3 per 1O 15 M . 
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Figure 2.6: Distribution of structures with respect to local density contrast, dn/dS(> 
1O 12 M ), for (3 = 2 and for Q m + Q A = 1 (left panel) and Einstein-deSitter (right 
panel) universes. Solid line: z = 0; dashed line: z — 1; dot-dashed line: z = 2; 
double-dot-dashed line: z = 3. The units of dn/dS are number of objects per Mpc 3 . 
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mass objects to be found in overdense regions was discussed by [13] and [12], who 
interpreted it by viewing halos today as progenitors of future larger-scale structures 
viewed at "high" or "low" redshift. 

In addition to the main peak at low \5\ , an additional, much lower and sharper peak 
can be seen right before the critical overdensity cutoff. This peak is the result of the 
change of the functional form of the conversion relation between linearly extrapolated 
and true density contrast close to virialization, when application of the spherical 
collapse model would lead 5 to diverge. The particular shape of the peak is an 
artifact of the recipe we adopted for dealing with the virialization regime, and carries 
no physical meaning (the shape of the peak is the shape of the high-5 end of dS/dS). 
However, since the boundary conditions we use for 5(5) and its derivative are physical, 
we do expect to have some form of local maximum at the high-5 end of the double 
distribution. Still, as discussed in the previous section, the effect of the details or even 
the existence of this local maximum on the physical quantities of interest is negligible. 

The high-<5 cutoff occurs at higher values of 5 in the concordance universe than 
in the Einstein-deSitter universe. This is a result of the different density contrast 
achieved at virialization by structures in the two different cosmologies. In the Einstein- 
deSitter case this density contrast is always 187T 2 , while in the concordance universe 
it is always higher and increases with time. At high redshifts, before the effect of A 
becomes significant, 5 C is very close to 187r 2 in the concordance universe as well, as 
can be seen in Fig. 2.3. 

Figures 2.4 and 2.5 show slices of the double distribution at various fixed values 
of 5, with z = and (3 = 2. Unlike the constant-mass slices, the constant-5 slices do 
not exhibit a peak (other than the global maximum imposed by the minimum-mass 
cutoff). Hence, there does not exist "most probable" mass at each given value of 
5. Hence, the hierarchical behavior of the number density of structures seen in the 
Press-Schechter mass function (which also exhibits no global maximum but instead 
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Figure 2.7: Distribution of density of matter inside collapsed structures with respect 
to local density contrast, dp/dd(> 1O 12 M ), for (3 = 2 and for f2 m + f2 A = 1 (left 
panel) and Einstein-deSitter (right panel) universes. Solid line: z = 0; dashed line: 
z — 1; dot-dashed line: z = 2; double-dot-dashed line: z — 3. The units of dp/d5 are 
M per Mpc 3 . 
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Figure 2.8: Fraction by number / n>un (solid line) and my mass f P;Un (dashed line) of 
objects of mass > 10 12 M living in underdense regions, as a function of redshift, for 
Q m + Qa = 1 (left panel) and Einstein-deSitter (right panel) universes. 
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diverges at low masses) also extends to the constant-5 slices of the double distribution. 
Figure 2.4 shows slices corresponding to relatively low values of \S\ (5 = —0.5,0,0.5 
and 3, close to the distribution peak in 5). At the high-mass end of the distribution, 
the abundance of objects increases with increasing 5, while in the low mass end of the 
distribution the trend is reversed, and the object abundance increases with decreasing 
5. This is in agreement with the behavior observed in the constant-m slices. 

Figure 2.5 shows slices corresponding to high values of 5 (5 = 10,20 and 30), 
farther from the distribution peak. In this case, the curves do not cross, and an 
increase of 5 simply results in an overall suppression of object abundance: structures 
of all masses are unlikely to be found overly clustered. This is because the final stages 
of collapse proceed rather quickly compared to the time spent around turnaround. 
The likelihood of a region observed in its late stages of collapse but before virialization 
is then low because the lifetime of this phase is small. 

In Fig. 2.6 we plot dn/d5(> 1O 12 M ) as a function of 1 + 5 for different values 
of redshift. It is striking that at z — 0, the distribution peaks at negative 5 values 
(around 5 = —0.6 in the concordance and —0.7 in the Einstein-deSitter universe), 
indicating that the most probable location for a collapsed object of mass > 10 12 M Q 
is an underdense environment. For the specific mass range, this trend is reversed 
by z — 1, when the preferred location of these objects is close to the universe mean 
(5 = 0). This time-evolution pattern is independent of cosmology, as it is present both 
in the concordance and the Einstein-deSitter universes, and appears rather to be a 
characteristic of the hierarchical nature of structure formation. Parameters of this 
distribution can be calculated using Eqs. (2.96) and (2.97) which, for the concordance 
cosmology and z = give (S) n = 0.43 and a$ jn = 4.36. The large value of the variance 
shows that the distribution is significantly broad. However, the positive value of the 
mean is an artifact of the asymmetric boundaries of the distribution and its long 
high-<5 tail. This is demonstrated by the notably different locations of the mean and 
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Figure 2.9: Slices of the double distribution function at m = 5.5 x 10 14 , M Q and for 
different values of the clustering scale parameter (3, for fl m + = 1 (left panel) and 
Einstein-deSitter (right panel) universes, plotted in linear scale. Solid line: (3 = 1.5; 
dashed line: /3 = 2; dot-dashed line: (3 = 3; double-dot-dashed line: (3 = 10. The 
units of the double distribution are number of objects per Mpc 3 per 10 15 M Q . 

the median. The value of the latter is 5 = —0.22, therefore more structures in this 
range reside inside underdensities. 

Figure 2.7 is the matter-density counterpart of Fig. 2.6, as it shows dp/d5(> 
1O 12 M ) as a function of 1 + 5 for the same values of redshift as in Figure 2.6. Again, 
at the current cosmic epoch, the distribution peaks at negative values of 5. Most of 
the virialized matter in the universe today appears to reside inside isolated objects 
rather than in clusters (note that decreasing the value of m min will only enhance 
this result since the trend towards isolation is more pronounced for the lower-mass 
objects). The trend of the peak with time (towards larger S for higher redshifts) is 
duplicated here as well. In particular, note that at present, a significant fraction of 
the mass lies in moderately underdense regions. Equations (2.98) and (2.99) give 
for this distribution (in the concordance cosmology and for z = 0), (8) p = 1.20 and 
0"<5,p — 6.23. The median of this distribution is at 5 = 0.20, a positive value. 

Finally, Fig. 2.8 shows the evolution with redshift of the fractions by number 
and by mass, / njUn and / ftm , of objects with m > 1O 12 M , living inside underdense 
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Figure 2.10: Distribution of structures of mass larger than 10 12 M with respect 
to local density contrast, dn/d5, for fl m + Q\ = 1 (left panel) and Einstein-deSitter 
(right panel) universes, at z = 0, and for (3 = 1.5(solid line), (3 = 2 (dashed line) and 
(3 = 10(dot-dashed line). The units of dn/dS are number of objects per Mpc 3 . 

regions. At high redshifts, when the mass of such objects is well above the exponential 
suppression cutoff, practically none of them are found inside underdensities. This 
trend is reversed as the redshift decreases. In the Q m + = 1 universe, an equal 
number of these structures are located inside underdensities by redshift 0.3 and by 
the current cosmic epoch, about 60% by number (but only 40% by mass) of these 
structures are located inside underdensities. 

Figure 2.9 demonstrates the effect of changing the clustering scale parameter on 
the double distribution. Slices of the double distribution along m = 5.5 x 1O 14 M 
are plotted (in linear function of S, and for f3 — 1.5 (solid line), 2 (dashed 

line), 3 (dot-dashed line) and 10 (double-dot-dashed line). The location of the peak 
appears to be extremely insensitive to the value of f3 for moderately low values. It 
very slowly moves towards 5 = with increasing /3, as it should (increasing (3 results 
in averaging the overdensity over increasingly large volumes). Note that even as (3 
approaches 1, the peak will not move towards 5 C , as a result of our correction for 
the central-object contamination. This makes our formalism particularly suitable to 
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study the properties of matter very close but outside a virialized structure (e.g. the 
local density of accreted gas). 

The effect of j3 on an integral quantity is shown in Fig. 2.10, which shows dn/d5(> 
1O 12 M ) for p = 1.5 (solid line), (3 = 2 (dashed line) and (3 = 10 (dot-dashed line). 
Again, the results are extremely insensitive to the value of j3, which gives us confidence 
about the robustness of the location of the peak of our distributions. 

As we have seen in this section, the correlations between the mass of collapsed 
structures and their environmental conditions are not only non-trivial, but also evolve 
dramatically with cosmic epoch (for a wide range of masses, the most probable 5 
changes sign between redshifts of zero and a few). It would therefore be desirable 
and of high interest to use the double distribution to investigate the effect of the 
environment on systems and phenomena the properties of which depend sensitively 
on the state of their surrounding space and matter. In the next chapter, we will use 
the environmental information encoded in the double distribution to assess the effect 
of the environmental factor on an important aspect of cosmic structure formation: 
cosmic accretion shocks. 
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Chapter 3 

Cosmic Accretion Shocks 



3.1 Overview 

In this chapter, we utilize the double distribution of dark matter halos which we 
derived and studied in the previous chapter, to study the effect of the environmental 
factor on the statistical properties of cosmic accretion shocks. In order to explicitly 
distinguish between effects of the environment and effects of the underlying distri- 
bution of accretor masses, we explore two variations of our model. In the first, all 
collapsed cosmic structures are assumed to live in a similar environment. The proper- 
ties of the accretion shock surrounding each such structure is then simply determined 
by the mass of the structure and the redshift of interest. The mass function of col- 
lapsed structures is, in this case, the Press-Schechter mass function. In the context 
of the second model, the distribution of collapsed structures with respect to both 
mass and environment is described by the double distribution. The properties of the 
accretion shock around each structure depends both on the mass of the structure 
as well as the environment of the structure (whether the structure resides inside an 
underdensity or an overdensity) . 

The formalism describing the properties of a single accretion shock around a cos- 
mic structure is presented in section 3.2. The Mach number of the accretion shock 
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surrounding each collapsed object is derived from the temperature jump across the 
shock. The temperature of the gas behind the shock is simply taken to be the virial 
temperature of the collapsed object. The temperature ahead of the shock depends on 
the model used each time. For the Press-Schechter-based model, it is simply the mean 
temperature of the intergalactic medium at the epoch under consideration 1 . For the 
double-distribution-based model, the temperature also depends on the density of the 
environment of the collapsed structure. The deviation from the mean intergalactic 
medium temperature is calculated by assuming that it is only due to adiabatic heat- 
ing (cooling) because of the relative local (de) compression with respect to the cosmic 
mean density. Once the Mach number of the accretion shock has bene determined, 
the accretion velocity (in the shock frame) is then calculated from it, and used to 
derive properties of the structure such as the accreted mass current and the kinetic 
energy crossing the shock surface per unit time. 

Combining the properties of a single shock with an underlying distribution of ac- 
cretors (either the Press-Schechter mass function or the double distribution of cosmic 
structures), in section 3.3 we derive the statistical properties of the population ac- 
cretion shocks. The calculated quantities include the distribution of accretor number 
density, accretor mass, shock surface, mass current, and kinetic power crossing the 
shocks, with respect to accretion shock Mach number. In addition, we calculate the 
cosmic history of the kinetic power and of the mass current integrated over shock 
Mach number, as well as the cumulative kinetic energy and number of baryons which 
have been processed through accretion shocks by each redshift. 

Our results are presented in section 3.4. We find that the shock environment 
alters the physical impact of shocks on the intergalactic medium, as well as the 

1 Although we include the effect of reionization (the cosmic epochs we consider are post- 
reionization and we take the temperature of the diffuse intergalactic gas at a density equal to the 
cosmic mean to be ~ 10 4 K), we do not include filament preheating in our models. Hence, we treat 
all accretion shocks as external shocks. 
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cosmic history of the shock population. The most prominent environmental effect 
is the development, at low redshifts, of a natural bi-modality in the distribution, 
with respect to shock Mach number, of the kinetic power crossing accretion shocks. 
Concerning the physical impact of shocks on the intergalactic medium, we find that 
the cumulative energy input of accretion shocks by redshift ~ 3 is comparable to 
the energy required to reionize the universe. In addition, more than a third of all 
baryons in the universe have been shocked in accretion processes by the present cosmic 
epoch. Finally, we comment on the components of the shock populations found in 
cosmological simulations [61; 63] with which our models are directly comparable. 

3.2 Properties of a single shock 

Throughout this chapter, we assume an adiabatic equation of state, and we con- 
sider all shocks to be non-radiative. We also assume that any individual collapsed 
object as well as its accretion shock are spherically symmetric. We will take the ac- 
cretion shock position around each structure to coincide with the virial radius of each 
structure. 

The Mach number of a shock, Ai, is defined as the ratio of the velocity of the 
accreted material in the shock frame to the adiabatic sound speed of the accreted 
material. The Mach number is related to the temperature jump across the shock 
through [74] 



where T± and are the pre-shock and post-shock temperatures correspondingly, and 
7 is the ratio of specific heats (assuming that this remains constant across the shock). 
For a 7 = 5/3 gas, Eq. (3.1) becomes 




(3.1) 



Tj. _ {5M 2 -1)(M 2 + 3) 
Tx ~ 16.M 2 



(3.2) 
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In the limit A4 S> 1 this equation is further simplified, 



The pre-shock temperature can be written in terms of the adiabatic sound speed of 
the pre-shock material c s i, 

m = ^4 , (3.4) 

7 

where \x is the mean molecular weight of the accreted gas, and m p is the proton mass. 
If we also take T 2 to be the virial temperature of the accreting structure, 

rp _r /xm p Gm 2 / 3 (4 7 r/ c p m , ) 1 / 3 (l + z) 

T 2 - T vir - — , (3.5) 

then the ratio T 2 jT\ becomes 

T 2 (m,z) r vir (m, z) 

( f \ 1/3 
= 2.7 x 10 3 fi m ( (1 + *) 

2 /3 /lr; 1 „-l\ 2 



m\ /15 kms 
x I — 



m 8 J V c si 



(3.6) 



where m is the object mass, f c = p v - ir / p m is the compression factor for a virialized 
object (which may vary with virialization redshift, depending on the cosmological 
model), z is the virialization redshift, h is the dimensionless Hubble parameter, and 
m§ = 5.96 x lO 14 /i -1 f2 m M is the mass included in a sphere of radius r§ = 8h~ x Mpc 
assuming the mean matter density inside the sphere to be equal to the cosmic mean. 

In Appendix D we compare this result with the Bertschinger similarity solution 
for an Q m = 1 universe [59] , and we find it to be in excellent agreement in the high-.M 
regime, where the Bertschinger solution is applicable. 

The surface area of a spherical shock around a structure of mass m is 

Si(m) = 47rr v (m, z) 2 (3.7) 
59 



where r v is the virial radius of the structure 



1.4/i _1 Mpc x 



(3.8) 



The mass current, defined as the rate at which mass crosses the surface of a single 
accretion shock around a structure of mass m at an epoch z, is 

Ji = ^ = 4vrr2(m)wip b (2;)(l + £ B ) 

= 47rr2(m)n b p c (2;)(l + 5 s )^c S)1 (3.9) 

where pb is the cosmic baryon density at the epoch of interest and (1 + 5 S ) is the 
density enhancement (with respect to the cosmic mean) just outside the shock. 

Finally, the kinetic power crossing the accretion shock around a single structure 
of mass m is 

„ dE 1 dm 9 1 T 9 , 
Pi = — = — — v? = - J x v\ . 3.10 
dt 2 dt 1 2 1 v ' 

Hence, 

P x = 27rr v (m) 2 fi b p c (^)(l + 5 s )M 3 c 3 sl . (3.11) 

3.3 Properties of the population of Cosmic Accre- 
tion Shocks 

The quantities we will use to describe the statistical properties of the population 
of cosmic accretion shocks are: 

• The "number distribution" of shocks with respect to Mach number. This is de- 
fined as the distribution of the comoving number density of accreting structures 
per logarithmic Mach number interval of their respective accretion shocks, 

dH M^-, (3.12) 



dliiM dM 
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with units number of structures per co moving Mpc 3 . 

The "surface distribution" with respect to Mach number. This is defined as the 
shock surface area per logarithmic Mach number interval per comoving volume 
under consideration, 

dS M^j (3.13) 



dlnM dM 

with units of Mpc -1 (since it represents a ratio of shock surface over space 
volume). 

The "accretor mass" distribution of shocks with respect to Mach number. This 
is defined as the distribution of the comoving mass density in accreting struc- 
tures per logarithmic Mach number interval of their respective accretion shocks, 

dp (3.14) 



dlnM dM 
with units of M per comoving Mpc 3 . 

The "mass current distribution" with respect to Mach number. This is defined 
as the comoving mass current density crossing shock surfaces of logarithmic 
Mach number between In M and In M + d In M , 

dJ M^- (3.15) 



dlnM dM 
with units of M Q yr _1 Mpc~ 3 . 

The "integrated mass current" , J, which is the comoving mass current density 
crossing shock surfaces of any Mach number at a given cosmic epoch, with units 
of M yr" 1 Mpc~ 3 . 

the "cumulative processed mass", f t Jdt, which is the total mass density pro- 
cessed by shocks of any Mach number since some initial epoch tj, expressed as 
a non-dimensional shocked baryon fraction. 



61 



• The "kinetic power distribution" with respect to Mach number. This is defined 
as the comoving kinetic power density crossing shock surfaces of logarithmic 
Mach number between In A4 and In M. + d In Ai , 

d\n.M dM 

with units erg s" 1 Mpc~ 3 . 

• The "integrated kinetic power" , P, which is the comoving kinetic power density 
crossing shock surfaces of any Mach number at a given cosmic epoch, with units 
of ergs -1 Mpc -3 . 

• The "cumulative processed kinetic energy" , Pdt, which is the total kinetic 
energy density processed by shocks of any Mach number since some initial 
cosmic epoch t- u with units of eV per baryon in the universe. 

In order to explicitly identify the environmental effects on the statistical properties 
of accretion shocks, we will use two different models to calculate these quantities. 
The first will assume that all structures are accreting gas of a single temperature, 
and that the population of accreting objects is well described by the Press-Schechter 
mass function. In this case, there exists a one-to-one correspondence between accretor 
mass and Mach number of the associated accretion shock. 

The second model assumes that the accreted gas has a distribution of densities and 
hence temperatures (where adiabatic heating and cooling are assumed to calculate 
the relation between local density and temperature). The distribution of accretors 
with respect to both mass and local over-(or under-) density is described by the double 
distribution of cosmic structures. In this case, there is a distribution of possible Mach 
numbers for the accretion shock around a structure of a given mass, depending on 
the local overdensity of the accreted material. 

Features of our results exclusive to the second model will then be an effect of the 
environment in which accreting structures reside. 
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3.3.1 Accreted Material of a Single Temperature 



If we assume that all collapsed objects accrete baryons of a single temperature, 
then at a given redshift, all objects with an accretion shock of a given Mach number 
will have the same mass, 



m = m(Ji4, z) 



4.2 x 10 9 /r 1 M m 



x 



18tt 2 



1/2 



QmfcJ ^15 kms -1 
(5M 2 -l)(M 2 + 3y 3/2 
16M 2 



-1 + z) 



-3/2 



(3.17) 



We will also assume that the mass distribution of collapsed objects can be described 
by the Press-Schechter mass function [1; 3] 

2 



dn 
dm 



m, z) 



2 Pm,0 



7r 



din a 



dlnm 



exp 



S c (z) 



2 [a(m)}' 



(3.18) 



where 8 c (z) is the linearly extrapolated overdensity of an object which collapses at 
redshift z, a(m) is the square root of the variance of the linearly extrapolated field 
smoothed at a mass scale m, and p m>0 is the cosmic mean matter density at the 
present time. 

In this case, the statistical quantities describing the population of accretion shocks 
can be derived in a straight forward way. The number distribution of shocks with 
respect to Mach number is 



dn(M, z) KA dn . dm . . 

- n - ] - r ^M- [m (M,z), z )—(M,,) 



where m = m(A4, z) is given by Eq. (3.17). 

Similarly, the surface distribution with respect to Mach number is 

dS 



(3.19) 



dlnM 



dn ^ 2 dm 
dm v dM 



(3.20) 
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while the accretor mass distribution is 



dp 



dn dm 
M m- — 



dlnM dm dM 

The mass current distribution in this model becomes 



dJ , . . , , , dn dm 
-(M, z) = M- 



(3-21) 



Ji 



(3.22) 



dlnM dm dM 

where J\(M,z) is given by Eq. (3.9) with 5 S = 3.13 (the overdensity factor in the 
Bertschinger solution (see Appendix D) just outside the shock), while the integrated 
mass current is 



J(z) 



, , dn dm 
M — - — 

M=1 dm dM. 



.MM 



(3.23) 



and the cumulative processed mass is 



t 



Jdt 



to ^ dn dm 
M — 7^ — 
dm dM 



Hi Jti JM=1 

which, in a concordance cosmology becomes 



J_ r Zi f°° dn dm 
H Jo Jm=i dm 9M 



JxdMdt 



(3.24) 



Jdt 



J\dz dM 



(i + z)y/n A + n m (i + z y 



(3.25) 



where H is the Hubble parameter and Z\ is the redshift corresponding to time t\. We 
use Zj = 10. 

Finally, the kinetic power distribution is 



dP , . . . , . dn dm 
:(M, z) = M 



Pi 



d\n.M dm dM 

where Pi(M,z) is given by Eq. (3.10), while the integrated kinetic power is 

dn dm 



(3.26) 



P(z) 



M 

M=1 dm dM 



P x dM 



(3.27) 
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and the cumulative processed kinetic energy is 

-to 



Pdt 



fj 



1 f Zi f°° dn dm 



Ho J Jm=i dm 



Pidz dM 



(3.28) 



3.3.2 Accreted Material of Varying Temperature 

In the second variation of our model, we wish to relax the assumption that the 
temperature of the accreted material is the same for all structures. Assuming that 
adiabatic heating or cooling is the only process that causes deviations of the temper- 
ature of the gas outside collapsed structures from its mean value, we can relate the 
local sound speed, c s , to the local overdensity or underdensity, 5, where 

S = Plocal ~ Pm . (3.29) 

Pm 

Since for adiabatic heating and cooling Cg oc p^cai' we get 

c s = c s , avg (5 + l) (7 - 1)/2 , (3.30) 

where c Sjavg is the "cosmic average" sound speed (the sound speed of the intergalactic 
medium at a density equal to the cosmic mean at the epoch of interest). This for- 
malism can accommodate the case where a process (such as reionization) heats the 
universe almost homogeneously, therefore increasing the average temperature (and 
consequently the cosmic average sound speed c Sjavg ). 

In order to make further progress and be able to calculate measures of the statis- 
tical properties of the population shocks in this approximation, we need an analytical 
model for the environment of collapsed structures. For this purpose, we will use 
the double distribution (DD) of collapsed structures with respect to mass and local 
overdensity which we derived and studied in the previous chapter. In the context of 
the double distribution, we defined the "local environment" of a collapsed structure 
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through the clustering scale parameter, (3. The clustering scale parameter is a free 
parameter in the DD model, and is defined so that the "environment" of an object of 
mass m be a surrounding region in space which encompasses mass (3m. The double 
distribution was found to be given by Eq. (2.22). 

In Eq. (2.22), Sg is the local linearly extrapolated overdensity (or underdensity) , 
which is related to the true (calculated from the spherical evolution model) overden- 
sity of the environment sphere including the local object, Si through the exact relations 
given in chapter 2, or through the useful approximation represented by Eq. (2.83), 
which is accurate at a better than 2% level throughout its domain for all cosmologies 
of interest. 

Because in the problem of cosmic accretion shocks we are interested in the prop- 
erties (density and sound speed) of the material right outside the shock surface, we 
will adopt a small value for the clustering scale parameter, (3 = 1.1. Our results, 
however, are not sensitive to the exact value of (3 since, as we saw in chapter 2, the 
properties of the double distribution (when calculated as a function of 5 rather than 
5e) depend very mildly on (3 for small values of j3. 

Distribution of Sound Speeds 

From the double distribution and Eq. (3.30), we can immediately calculate the 
number density of collapsed objects of mass m > m min embedded in a medium of 
local sound speed between c s and c s + dc s , 



Similarly, the density of matter in collapsed objects of mass > m m j n embedded in a 
medium of local sound speed between c s and c s + dc s is 




(3.31) 




(3.32) 
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In Eqs. (3.31) and (3.32), d5/dc s can be found from Eq. (3.30), 

(3.33) 



d8 2 1 / c, \^ 



dc s 7 1 C Sjav g \C sav g^ 

Hence, if we use C s to denote the sound speed in units of the cosmic average sound 
speed, 

C s = -^, (3.34) 

Cs,avg 



dn 2 2=2 r°° dn 

Cr 1 / -r^dm, (3.35) 



Eq. (3.31) becomes 

dC s 7 - 1 ~ ° J m=mmin dmd5 

and Eq. (3.32) becomes 

dp 2 „S=% f 00 dn 



Shock Properties 



-C 1 ^ 1 / m dm (3 36) 

dC s 7 - 1 7m=m min dmcM 



In this second variation of our model, the mass of a collapsed object with an 
associated accretion shock of Mach M. is also dependent on the local overdensity S, 



m = m(M,6,z) = 4.2 x lO 9 /*" 1 M 



187T 2 N 1/2 



in 



L?' avg 1 ) 3 (^+n 3(7 ~ 1)/2 

\15 kms V 

(5.M 2 -1)(.M 2 + 3) 1 3/2 /2 



Note that, when 5 = 0, this equation is reduced to Eq. (3.17), as it should. 
Hence, the number distribution of shocks with respect to Mach number is 

dn . , f 5c , „ dn r ,.,„,,. , 9m 



(3.37) 



(3.38) 



where m(M, 5, z) is given by Eq. (3.37). 
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Similarly, the surface distribution with respect to Mach number is 
dS 



d^M = M Lj^dmTs dM 



(3.39) 



while the accretor mass distribution is 

dp 



dlnM 

The mass current distribution is 

dJ 



, Sc , . dn dm 
Ml dom- 

5=-l 



dmdS dM 



(3.40) 



8,z 



/ SC JIT dU dm 

M I d5J r 

8=-l 



(3.41) 



d\nM Js=-i dmdS dM 

where Ji(M, 5, z) given by Eq. (3.9). In this case, 5 S = 5, the local overdensity outside 
the shock as given by the DD. The integrated mass current is then 

dn dm 



J(z) 



oo rS c 

/ dMdSJi 
M=iJs=-i dmdd dM 



and the cumulative processed mass is 



t 



Jdt 



o r oo f s c dMdSdzJ^^l 



z t J M=l J -I 



h (i + z )y/n A + n m (i + zy 



Finally, the kinetic power distribution is 
dP 



d\nM 



5c dn dm 

M I d5Py 
S=-l 



dmdS dM 



with Pi(M,S, z) given by Eq. (3.10), while the integrated kinetic power is 
P(z) = 



/ dMddPy 

M=l JS=~1 



dmdb dM 



and the cumulative processed kinetic energy is 



o /-oo f 5 c dMdbdzPx 



to 



Pdt 



zi J M 



't; 

dn dm 



dmdS dM 1 8,z 



i J -i H (l + z) v /Q A + Q m {l + 



(3.42) 



(3.43) 



(3.44) 



(3.45) 



(3.46) 
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Figure 3.1: Distribution of number of objects (left panel, in units of objects per 
comoving Mpc 3 ) and mass density (right panel, in units of M Mpc -3 ) per C s interval, 
for objects with mass > 10 12 M , in a WMAP concordance universe. Solid line: z = 0; 
dashed line: z — 1; dot-dashed line: z = 2. 

3.4 Results 

3.4.1 Distribution of Environmental Sound Speeds 

In Fig. 3.1 we plot the distribution of sound speeds of the material in which 
collapsed objects of mass greater than 10 12 M Q are embedded, for different redshifts. 
The left panel shows the comoving number density of objects per interval of the 
dimensionless sound speed C s , dn/dC s , and the right panel shows the mass-density 
counterpart dp/dC s . All curves correspond to a concordance Wilkinson Microwave 
Amsotropy Probe (WMAP) universe (<r 8 = 0.84, h = 0.71, fi m = 0.27 and fi b = 0.04, 
[24]). As expected from the results of chapter 2, which show that the most probable 
density contrast of the material surrounding collapsed structures decreases with time, 
the most probable sound speed of the material surrounding collapsed structures of 
mass > 1O 12 M increases with redshift. The overall suppression of the curves with 
increasing redshift is due to the fact that fewer structures more massive than a fixed 
cutoff have had enough time to collapse at large redshifts. At the present epoch, the 
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Figure 3.2: Number distribution of shocks per logarithmic Mach number interval for 
the Press-Schechter-based model (left panel) and the double-distribution-based mod- 
els (middle and left panels). The middle panel models implement a halo mass cutoff 
identical to that of the Press-Schechter case, while the right panel model assumes no 
mass cutoff. The units of the vertical axes are number of objects per comoving Mpc 3 . 
Solid line: z = 0; dashed line: z = 1; dot-dashed line: z = 2. 

most probable sound speed for the masses under consideration is only a fraction of the 
cosmic average sound speed, due to the fact that structures in this range of masses 
are more likely to be found inside underdensities. Note however that the picture 
presented here is incomplete, as shock heating in filaments will further modify the 
sound speed distribution, favoring higher values of the sound speed. 

3.4.2 Properties of Accretion Shocks: Effects of the Local 
Environment 

In this section we again assume a WMAP concordance universe. We focus on 
post-reionization redshifts, hence we assume a cosmic average sound speed of 15 km/s, 
corresponding to a temperature of ~ 10 4 K, for a fully ionized plasma with \i = 0.59 
(25% He by mass). 

In Fig. 3.2 we plot the number distribution of shocks with respect to Mach number, 
dn/dlnA4, for the Press-Schechter-based (hereafter PS) model (left panel) and two 
implementations of the double-distribution-based (hereafter DD) model (middle and 
right panel). 
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The distribution in the PS model monotonically increases for decreasing M. since 
as we have discussed in this case there is a one-to-one correspondence between mass 
and Mach number, hence the number distribution of shocks simply follows the Press- 
Schechter mass function. The requirement M. > 1 imposes then a cutoff in the mass 
of the dark mater halos that can act as a host to an accretion shocks. This mass cutoff 
is the mass of the object whose virial temperature exactly equals the temperature of 
the diffuse gas, and it is m min = 8.1 x 1O 9 M at z = 0, m min = 3.8 x 1O 9 M at z = 1 
and m m in = 2.2 x 1O 9 M at z = 2. The decrease of m m m with increasing redshift is 
due to the fact that the virial temperature of an object of a fixed mass is larger if the 
object virializes at high redshift than if it virialized today. In an Einstein-deSitter 
Universe this is a simple result of the higher cosmic matter density in the past (the 
compression factor of a virialized object is constant in time, so since the cosmic 
density is larger in the past, the virial density which is the compression factor times 
the cosmic density is also larger in the past, resulting to a higher virial temperature). 
In the concordance universe, there are two competing factors mediating this effect. 
The cosmic matter density is still larger in the past, however the compression factor 
is now itself redshift- dependent, and increases with decreasing redshift. However, 
this increase of the compression factor is not steep enough to counteract the decrease 
in cosmic density, and still the virial density decreases with decreasing virialization 
redshift . 

The middle panel of Fig. 3.2 shows a variation of our DD model which includes 
only objects of a mass larger than the mass cutoff of the PS model in the correspond- 
ing redshift. The mass cutoff, in combination with the double-distribution-imposed 
distribution of pre-shock sound speeds for each accretor mass, results in a M. number 
distribution which peaks at Ai < 4, a position which is defined by the combination 
of the virial temperature at the mass cutoff (which is the most populated available 
mass bin) and the sound speed corresponding to the most probable environment (as 
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given by the double distribution) at that particular mass. If the additional constraint 
M. > 1 did not exist, then the area below the curves corresponding to the same 
redshift in the left and middle panels would be the same, as it would represent simply 
the total number of objects with mass larger than the mass cutoff. However, some of 
these objects are "lost" in the double distribution model because they are embedded 
in material of temperature equal or larger than their virial temperature and hence 
they cannot harbor shocks. This effect is more pronounced at higher redshifts, while 
at z = where very few objects are "lost" from the distribution as described above, 
the areas below the two solid lines in the left and middle panels are almost equal. 

The right panel of Fig. 3.2 shows the DD model where the only restriction in the 
participating objects is the requirement that M. > 1. In this case, the distribution 
does not turn over, and it is in fact at an overall higher level that the corresponding 
PS model. This is because the DD model tends to move smaller objects to higher M. 
bins, because in the double-distribution picture smaller objects preferentially reside 
inside underdense regions. Therefore, smaller objects accrete cooler gas and harbor 
shocks of higher M, than in the Press-Schechter case. When no mass cutoff is im- 
posed, smaller objects (that are excluded in the Press-Schechter picture because their 
virial temperatures are below the cosmic mean temperature) "leak" inside the M. 
distribution, because in the double-distribution picture they can lie inside an under- 
dense region (void) and accrete material cooler than the cosmic mean and cooler than 
their own confined gas. As a result, a larger total number of objects participate in 
the double distribution and the overall level of the curves is higher than in the PS 
model. 

Note that the effect described above is moderated by the fact that some objects 
also "leak" out of the distribution, because they lie in overdense regions where the 
diffuse gas has a temperature larger than their virial temperature (it is this effect 
that caused the "loss of signal" that is observed in the middle-panel curves). 
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In all panels of Fig. 3.2 the amplitude of the distribution increases with increasing 
redshift. This is because the number distribution plotted here is dominated by the 
low-mass objects, which are merged into larger halos as time progresses. In Fig. 3.3 
we plot the double distribution without cutoff model (right panel of Fig. 3.2) using 
a logarithmic scale for the vertical axis, so that the behavior of the large-mass (and 
consequently high-A^) objects can be better illustrated. At these high A4s, the 
behavior of the distribution with redshift is reversed, as more high-mass objects are 
formed as time progresses. 

Figure 3.4 shows the distribution of spatial density of shock surface per logarithmic 
Mach number interval dS/d In M. , with (dS '/ 'din M)d\\i M. being the total surface of 
shocks in a certain comoving volume V with logarithmic Mach number between In Ai 
and In A4 + d\nA4 divided by that volume. 

This quantity is frequently used to characterize the statistical properties of shocks 
in cosmological simulations as it can be calculated without need for identification of 
collapsed structures and assignment of shock-hosting gridpoints to specific structures. 
The left panel shows this distribution for the PS model while the middle and right 
panels correspond to the DD models, with (middle panel) and without (right panel) 
mass cutoffs. The behavior of this quantity is similar to the number distribution with 
respect to M. in that it is dominated by the low-mass structures, and therefore the 
existence of a low-mass cutoff significantly affects its qualitative behavior. 

Although the number and surface distributions of shocks with respect to Mach 
number are dominated by the low-mass objects and their corresponding low-Mach- 
number shocks, the accretor mass, the mass accretion rate and the energetics of shocks 
are dominated by high-mass objects. This can be immediately verified by simple 
analytic arguments. We consider for simplicity the Press-Schechter-based variation 
of our model. For masses high enough that the primordial density fluctuation power 
spectrum can be regarded as a power law but low enough that the exponential mass 
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Figure 3.3: Distribution of number density of objects per logarithmic Mach number 
interval for the double-distribution-based model without a mass cutoff, plotted in 
logarithmic scale. The units of the vertical axes are number of objects per comoving 
Mpc 3 . Solid line: z = 0; dashed line: z — 1; dot-dashed line: z — 2. 
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Figure 3.4: Distribution of spatial density of shock surface per logarithmic Mach 
number interval (dS/s In A^)for the Press-Schechter-based (left panel) and the 
double-distribution-based (right panel) models. The units of the vertical axes are 
comoving Mpc -1 . Solid line: z = 0; dashed line: z — 1; dot-dashed line: z — 2. 

cutoff is not affecting the results, dn/dm oc m~ 2 . In addition, M oc m 1 / 3 (in the 
high-AI limit), while r v oc m 1//3 . Hence, J\ oc m and Pi oc m 5//3 . Equations (3.19)- 
(3.26) then give for the low-mass dependence of the various statistical properties of 
the shock population, 



dn 



dlnM 

dS 

dlnM 

dp 

dlnM 
dJ 

dlnM 
dP 



oc m 



oc m 



-1/3 



oc m 



o 



oc m 



o 



oc m 



2/3 



(3.47) 



dlnM 

accounting for the difference in the low-mass behavior of different distributions, in par- 
ticular the divergent behavior of the number distribution and the convergent nature 
of the mass, mass current and kinetic power distributions. Note that the low-mass 
suppression is somewhat stronger than what is predicted by the simple arguments 
above, due to the deviation of M(m) from m 1//3 for M — > 1. 

In Fig. 3.5 we plot the accretor mass distribution with respect to Mach number, 
i.e. the mass in objects with accretion shocks of logarithmic Mach number between 
log.M and log.M +dlogM versus M. All models exhibit similar behavior, with the 
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Figure 3.5: Accretor mass distribution (mass density in accretors per logarith- 
mic Mach number interval of associated accretion shock, dp/ d\aM)ioi the Press- 
Schechter-based (left panel) and the double-distribution-based models, with a Press- 
Schechter like mass cutoff (middle panel) and without an explicit mass cutoff (right 
panel) . The units of the vertical axes are M co moving Mpc" 3 . Solid line: z = 0; 
dashed line: z — 1; dot-dashed line: z = 2. 

high-AI objects dominating the distribution. In this case, the presence or not of a 
mass cutoff in the DD models does not have a significant effect, as it only affects the 
low-.M end of the distribution. The location of the peak is at similar values of M. in 
all models, and is an effect of the high-mass exponential cutoff in collapsed accreting 
structures. The effect of a varying environmental overdensity and sound speed in this 
case is to increase the height as well as steepness of decline toward both higher and 
lower Mach numbers of the distribution peak. This is because the double distribution 
tends to move lower-mass objects (which generally reside in underdensities) towards 
higher A4s, and higher-mass objects (which reside in overdensities) towards lower 
.Ms, causing a greater concentration of accretor mass close to the peak. 

In Fig. 3.6 we plot the mass current distribution, dJ/d In Ai. Again, the left panel 
shows the results for the PS model, the middle panel represents the DD model with 
a Press-Schechter mass cutoff and the right panel is the DD model without a mass 
cutoff and with the only restriction being M, > 1. 

The environmental effects are more pronounced in this case, while the existence or 
not of a mass cutoff plays no significant role in the properties of the DD model, again 
due to the overwhelming dominance of high-mass objects in the total mass processing 
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Figure 3.6: Mass current distribution (spatial density of mass current per logarith- 
mic Mach number interval (dJ/d In .M)for the Press-Schechter-based (left panel) and 
the double-distribution-based (right panel) models. The units of the vertical axes are 
M comoving Mpc~ 3 yr. Solid line: z = 0; dashed line: z — 1; dot-dashed line: z — 2. 

rate through accretion shocks: the existence or not of a low-mass cutoff cannot affect 
the mass current distribution because the \ow-A4 objects process only a very small 
fraction of the total accreted mass. 

The effect of taking into account the environmental overdensity distribution using 
the DD in this case is two-fold: On the one hand, it spreads out every mass bin to a 
larger Ai range as in the previous cases. On the other hand, it also adjusts the local 
density of the material just outside the object, which in turn affects the local value 
of the mass current. 

The location of the mass current distribution peak and its evolution with red- 
shift do not change appreciably between the two models. The location of the peak 
represents the high-mass exponential cutoff in the Press-Schechter mass function (the 
largest mass scale which has collapsed by a certain epoch), modulated by the nonlinear 
mass-.M relation. However, the DD model declines more sharply towards lower A4. 
This is because in the double-distribution picture, low-mass, low-.M objects reside in- 
side underdensities and hence they can process less mass than in the Press-Schechter 
case. 

Finally, the suppression of the distribution amplitude with decreasing redshift in 
all models is a result of the reduction of the mean cosmic density of the accreted 
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Figure 3.7: Kinetic power distribution (spatial density of kinetic power processed 
by accretion shocks per logarithmic Mach number interval, dP/dlnAi) for the Press- 
Schechter-based (left panel) and the double-distribution-based (right panel) models. 
The units of the vertical axes are 10 40 ergs s -1 comoving Mpc~ 3 (left axis) and 10 51 
ergs per century per comoving Mpc 3 (right axis). Solid line: z = 0; dashed line: 
z = 1; dot-dashed line: z = 2. 

material. The result is more pronounced in the DD models since, in addition to the 
reduction of the mean density, there is also a shift of the most probable environmental 
density for structures of a given mass towards regions increasingly underdense with 
respect to the cosmic mean as time increases. 

In Fig. 3.7 we plot the kinetic power distribution with respect to Mach number, 
for the PS model, and the two variations of the DD model. As in the case of the mass 
current distribution, the kinetic power distribution is also dominated by objects of 
high mass and Mach number, and the presence of a lower-mass cutoff does not have 
a significant effect on the DD models. 

A most striking environmental effect in this case is that, as redshift decreases, 
a second peak separates out in the DD model, which is absent in the PS models. 
The presence of this second, high-.M peak is due to the gradual shift of increasingly 
massive structures towards underdense environments. To better demonstrate this 
effect, we plot in Fig. 3.8 the components of the distribution produced by structures 
in significantly overdense environments (dot-dashed line, 5 > 5 tau , where <5 ta is the 
overdensity of a perturbation turning around at the specific epoch), and structures in 
environments of 5 < <5 ta (dashed line), for the DD model without a mass cutoff and 
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Figure 3.8: Distribution of spatial density of kinetic power per logarithmic Mach 
number interval (dP/s In M) the double-distribution-based model without a mass 
cutoff. Left panel: z = 2; middle panel: z — 1; right panel: z = 0. Solid line: 
overall distribution; dot-dashed line: contribution from structures with environmental 
overdensities between £ ta and S v (between turnaround and virialization overdensities); 
dashed line: contribution from structures with environmental overdensities 5 < 5 ta . 

for three different redshifts. The contribution from the low-density environments is 
distinct in Ai space and increasing with decreasing redshift. 

Figure 3.9 shows the evolution of the integrated kinetic power over shocks of any 
Mach number, P, for redshifts between 10 and 0. The solid line corresponds to the 
DD model and the dashed line to the Press-Schechter models. In the left panel, P{z) 
is plotted as a function of z, and is seen to peak near z ~ 1 for both models, when the 
increase with time of the number and mass of collapsed structures is balanced by the 
decreasing mean density of the accreted material due to cosmic expansion. The effect 
is more pronounced and the peak occurs at a slightly higher redshift in the DD model 
where environmental effects are accounted for, since the decrease in the cosmic mean 
density is accompanied by a decrease in the most probable overdensity with respect 
to cosmic mean. The overall level of the curve in the DD model is also a factor of 
~ 2 higher compared to the PS model at redshifts > 1, because environmental effects 
further enhance the contribution of the larger structures which are already favored 
due to their size and accelerating potential. However, at lower redshifts, when in the 
DD model the contribution of the smaller structures becomes comparable to that of 
the larger structures, the integrated kinetic power predicted by the two models tends 
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Figure 3.9: Integrated kinetic power over shocks of any Mach number, P, for the 
Double-Distribution (solid line) and the Press-Schechter (dashed line) models. Left 
panel: redshift history of P for z < 10. Right panel: J Pdt in units of eV per baryon. 
The horizontal line in the right-panel plot corresponds to 13.6 eV per baryon. 

to converge. 

The right panel shows the redshift history of J Pdt (the cumulative processed 
kinetic energy) in units of eV per baryon in the universe (as opposed to per shocked 
baryon). Again, the solid line is the DD model while the dashed line is the PS model. 
The horizontal line in this plot corresponds to 13.6 eV per baryon. From the location 
of the intersection of the horizontal line with the J Pdt curve, we can conclude that 
the energy processed by accretion shocks alone by redshift z ~ 3 (~ 2 in the PS 
model) is of order of magnitude comparable to the energy required to reionize the 
universe even in absence of other sources of energy 2 . 

Finally, in Fig. 3.10 we plot the integrated mass current over shocks of any Mach 

number, J, for the DD (solid line) and the PS (dashed line) models. The left panel 

shows the redshift history of J for z < 10, while the right panel shows the redshift his- 
2 Note that this order-of-magnitudc argument is meant to give a feeling about the amount of 
energy processed by shocks as compared to other energy inputs in the IGM. If one wanted to consider 
shocks as an actual reionization mechanism, a detailed modelling of the reionization process would 
be required. 
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Figure 3.10: Integrated mass current over shocks of any Mach number, J, for the 
Double-Distribution (solid line) and the Press-Schechter (dashed line) models. Left 
panel: redshift history of J for z < 10. Right panel: J Jdt expressed as the fraction 
of baryons in the universe which have been processed by accretion shocks. 

tory of J Jdt, the cumulative shocked gas mass, expressed as the fraction of baryons 
in the universe which have been processed by accretion shocks. By the current cos- 
mic epoch, a fraction between 40 — 50% of baryons have already been processed by 
accretion shocks, while ~ 10% of the baryons have already been processed by shocks 
by z ~ 3. 

3.4.3 Comparison with Cosmological Simulations 

Caution should be exercized when comparing these calculations to results of cos- 
mological simulations studying the properties of cosmic shocks, to ensure that corre- 
sponding quantities are being compared. One should keep in mind that the shocks 
studied in this chapter (a) are accretion shocks only, while effects of merger and fil- 
ament shocks and their contribution to shock statistics and energetics have not been 
included and (b) are all considered to be external shocks, accreting material which 
has never been shocked before. 

In [63], the cosmological simulation results are presented by classifying shocks 



81 



as external or internal, depending on the temperature of the accreted material. In 
both the kinetic and surface distributions with respect to Mach number, used by 
the authors to describe the statistical properties of cosmic shocks, external shocks 
are dominated by filament shocks, while internal shocks are dominated by merger 
shocks. The population of accretion shocks is presumably divided between these 
two categories, depending on whether structures accrete void or filament material. 
Accretion shocks are expected to lie in the high Mach number end of both internal 
and external shock distributions. Hence, no direct comparison between our model and 
these results is possible without the inclusion of a detailed model of filament- shocked 
material and the way it modifies the properties of the accretion shock population. 
However, properly adjusting our model to the [63] simulation parameters, we find 
that if we consider the entire population of accretion shocks to be either external (as 
in the results of the previous section) or internal with mean density and temperature 
values for the filament gas taken from the [63] results, our curves lie, as expected, 
below their overall curves for external or internal shocks, respectively, for all relevant 
quantities, although they constitute an appreciable fraction of them. 

In [61], the simulated shocks are instead divided in accretion and merger shocks, 
the accretion shocks being assigned to specific structures. The principal quantity 
used to statistically describe the shock population itself is the number distribution of 
shocks with respect to Mach number. However, a direct comparison is not appropriate 
in this case either, because part of the simulated population of accretion shocks resides 
inside filaments, which introduces a bi-modality in the number distribution, especially 
pronounced in the Q m = 1 simulation. In addition, the properties of the shocks even 
around structures which accrete principally pristine material are modified by lateral 
accretion of hotter gas from the filaments. Finally, the combination of selection 
criteria for the collapsed structures inhibits the identification of a sharp mass cutoff 
to use in our analytical models. 
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Overall, although the models we have used here can reveal the effect of a certain 
class of environmental factors (local overdensity and associated change in tempera- 
ture) on the statistics of cosmic accretion shocks, it becomes clear that for a detailed 
comparison with simulations and observations, a treatment of the second important 
environmental factor, i.e. the filamentary structure of the universe, is needed. 

A result of this work with direct relevance to the interpretation of cosmological 
simulations studying the properties of cosmic shocks is the effect of the existence of 
an explicit mass cutoff to the properties of different shock distributions, since mass 
cutoff's are always present in simulations as a result of finite mass resolution. We have 
seen that the number distribution as well as surface distribution of shocks with respect 
to Mach number are dominated by the contribution of low-mass objects and hence 
they are appreciably affected by the presence of a mass cutoff. The rest of our shock 
distributions however, such as the mass current distribution and the kinetic power 
distribution are dominated by high-mass objects and are unaffected by the presence 
of a mass cutoff. This is consistent with the findings of [63], who performed a series 
of simulations with increasing mass and spatial resolution, and found convergence 
for their results on kinetic power distribution but not for the surface distribution. 
However, the physical impact of shocks to their environment is better represented by 
the latter class of distributions, which are unaffected by the value of the mass cutoff, 
and hence the lack of convergence in the number-dominated distributions need not 
decrease confidence in the relevance and robustness of the physical output of such 
cosmological simulations. 
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Chapter 4 



Discussion 



In this thesis, we have derived and explored a new tool for the analytical study 
of hierarchical structure formation: a double distribution of the number density of 
cosmic collapsed structures with respect to mass and local overdensity. We have 
done so by introducing a clustering scale parameter /3 > 1, which we use to associate 
with each collapsed object of mass m a larger environment of mass (3m. The scale 
parameter (3 can be expressed as a function of the number of virial radii included 
in the local environment of each structure. We found that for reasonable values 
(3 ~ 2, the shape of the distribution does not depend sensitively on this parameter. 
Integration over linearly extrapolated overdensity returns the original Press-Schechter 
mass function, independently of the value of (3. 

We have presented the double distribution in terms of the true, physical, nonlinear 
density contrast 5. However, in calculating the distribution it is useful to identify 
regions using instead the overdensity obtained via linear analysis, 5, extrapolated to 
the present epoch. A useful fitting function was given for the 5 — 5 conversion. 

The double distribution is useful because it allows us to have an explicit analytical 
if approximate description of the environment in which collapsed objects of all masses 
reside. Using the tools we have developed, it can be readily calculated for any flat 
cosmology, and evaluated at any epoch. Consequently, it offers new insight into the 
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growth of structure as well as the present distribution of collapsed objects. 

We have evaluated the double distribution and some of its integral moments for 
both a concordance cosmology and an Einstein-de Sitter universe. Some key results 
are that at any redshift, the double distribution is dominated by a peak which shifts 
in mass but is always at a relatively low value of \8\. For each mass, there is a most 
probable 5, which increases with structure mass. Moreover, at the present epoch in 
the concordance universe, the most probable environment is a modest ■unc/eniensity, 
for all objects below about 10 14 M Q ; thus, underdensities are preferentially populated 
by low-mass objects. Finally, the fraction of mass in underdensities increases with 
time, and in the concordance cosmology the present underdense mass fraction in 
objects of M > 10 12 M is about 40%. These trends can be understood in terms 
of hierarchical clustering in which overdense regions are the site of vigorous merging 
that clears out low-mass objects, which then find their last refuge in voids. 

In addition, we have utilized the double distribution of cosmic structures to in- 
vestigate analytically the effect of environmental factors on the properties of cosmic 
accretion shocks around collapsed structures. For this purpose, we have explored two 
different models for the cosmic shock population. The first used the Press-Schechter 
mass function to describe the underlying population of collapsed, accreting objects. 
All such objects were assumed to accrete material of the same density and tempera- 
ture. This was our "control" model, which did not include any environmental effects. 
The second model used the double distribution of collapsed structures [75] to describe 
the distribution of accreting objects with respect to both their mass and local envi- 
ronment over density or underdensity. The overall mass distribution of objects is the 
same as in the first model, as the double distribution integrates to the Press-Schechter 
mass function. 

We found that the number and surface distributions of shocks with respect to Mach 
number are dominated by the contribution from low-mass objects and hence peak at 
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low Mach numbers for both models. The contrary is true for distributions describing 
the physical impact of accretion shocks on their environment, such as the mass current 
and kinetic power distributions, which are dominated by the properties of primarily 
high-mass objects, and hence peak at high Mach numbers. The distribution peaks 
are more pronounced when environmental effects are taken into account, as they tend 
to move objects of both high and low masses towards the peak. 

Perhaps the most striking effect of accounting for the environmental factor is the 
separation of a second, high Mach number peak in the kinetic power distribution at 
low redshifts, due to an increasing number of higher-mass structures concentrating 
inside underdensities. This double-peaked behavior is present despite the fact that 
the effect of filament heating of the accreted gas was not included in our calculation. 
The latter process is expected to further complicate the features of the kinetic power 
distribution. 

The integrated kinetic power processed by shocks peaks at a redshift of ~ 1, and 
the effect of the local environment is to increase the overall level of the processed 
energy at high redshifts by a factor of ~ 2 as well as to move the peak of the kinetic 
power history towards slightly higher redshifts. We found that accretion shocks alone 
have processed by z ~ 3 energy comparable to that required to reionize the universe. 

The integrated mass current history (i.e., the net baryonic mass processed through 
shocks) peaks at earlier epochs compared to the kinetic power redshift, with the peak 
in both models occurring at z ~ 3. By the current epoch, the baryon fraction shocked 
in accretion shocks alone is between 40 — 50%. Since this material represents baryons 
which can condense to form galaxies and stars, this fraction represents an upper limit 
to the baryon fraction in galaxies. 

The population study of cosmic accretion shocks presented in the second part 
of this thesis is only one of many possible applications of the double distribution of 
cosmic structures. With its capability to treat both underdensities and overdensities, 
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and collapsed as well as diffuse regions, the double distribution can be used to build 
analytical models for many problems which are currently attracting much interest 
in the context of cosmological structure formation. The filamentary structure of the 
universe and the dependence of merger histories on environmental factors are only 
two such problems that we are currently pursuing using the tools developed in this 
thesis. 
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Appendix A 



Limits of the Double Distribution 



In this appendix, we examine the behavior of the double distribution in the limiting 

cases (3 — > oo and j3 — > 1. 



A.l Behavior of the Double Distribution in the 



limit (3 — ► oo 



In order to find the behavior the double distribution as (3 — > oo, we recall that, 
because S(m) decreases monotonically with m, its limit in the infinite (3 regime will 
be 

lim S((3m) = . (A.l) 
Then, using the notation of the previous section, 



lim — (m, 5g, (3, a) = lim —(Si,S 2 ,Se,a,m). 
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The limit of a unit-area Gaussian when its width vanishes is the Dirac delta-function 
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Using this result, we get 
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However, the 5g— domain of the double distribution is between — oo and <5o jC , and 
therefore the value Sg = 2<5 , c is outside its domain. Hence the second Dirac delta- 
function is always zero, and 
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proportional, as expected, to a Dirac delta-function centered at Si = 0. 



A. 2 Behavior of the Double Distribution in the 



limit j3 



Denoting S(j3m) by 5*2 and letting <fi = S(m)/S(/3m), we seek the behavior of the 
double distribution in the limit — > 1 or — > 1. Defining 
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we can write 
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Appendix B 



Derivation of the Press-Schechter 



Mass Function Prom the Double 



Distribution 



Using 5*i to denote S(m) and S 2 for S((3m) we have: 
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Then, performing the transformation 8' — > —5' in the second integral, we get 
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The final result is the Press-Schechter mass function formula, independently of the 
value of (3. 
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Appendix C 



Vacuum Integrals 



C.l The incomplete vacuum integral of the first 
kind Vi 

C.l.l Definition 



with domain < r < 1 and jj > 2. 
C.1.2 Properties 

Physically, V\ (r, /x) is proportional to the time required by a perturbation of normal- 
ized curvature parameter k/oj 1 ^ = (fi + l)//i 2 / 3 to achieve a size a p = rap^n/uj 1 / 3 ) 
before turnaround. Its asymptotic behavior for r 1 is 



We define the incomplete vacuum integral of the first kind as 
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while for /x ^> 1 it is 
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In the case of the Eddington perturbation (/i = 2), we can derive a closed-form 
expression for Vi: 



Vi(r,2) 
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When r = 1, the value of Vi(l,/i) is the complete vacuum integral of the first kind, 
which is a function of \i alone. Physically, the complete vacuum integral of the first 
kind is proportional to the time required for a perturbation of curvature parametrized 
by /i to reach turnaround. The derivative of Vi(l, /x) appears in the calculation of the 
derivative d5 Q /d5, in the 3rd line of Table 2.4, and it is 



djj, 



y/1 — x{—X 2 — X + /i) 3 / 2 



(C.5) 



C.2 The hyperbolic vacuum integral of the first 
kind Tii 



C.2.1 Definition 

We define the hyperbolic vacuum integral of the first kind as 

3 f *\J~xdx 



2 J v /(l + x)(x 2 -a; + w) 
with domain < r < 00 and w > 1/4. 
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C.2. 2 Properties 

Physically, 7i\(r,m) = is proportional to the time required by a perturbation of 
normalized curvature parameter k/lu 1 / 3 = (1 — zo)/vj 2 / 3 to achieve a size a p = 
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rdp^n/uj 1 / 3 ). Its asymptotic behavior for r <C 1 is 

r<l 1 



Wifct^'w -^r 3 / 2 (C.7) 



while for r ^> 1 it is 

H a (r, n?) ^ C(o7) + - In far 2 -r + w + 2r - 1) (C.8) 

where C(w) is a function dependent only on w. In the case of a flat [w = 1) 
perturbation, Hi can be integrated immediately to give 
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C.3 The incomplete vacuum integral of the second 
kind V2 

C.3.1 Definition 

We define the incomplete vacuum integral of the second kind as 

V2(r '^ = 4 7 (l- X )3/2 ( _ a; 2_ x + ^ ) 3/2. ( C - 10 ) 

with domain same as for Vi(r, /i). 
C.3. 2 Properties 

The incomplete vacuum integral of the second kind is related to Vi(r, /x) through 

-V 1 (r, / u)= / u 2 / 3 V 2 (r, / u), (C.ll) 
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In the case of the Eddington perturbation (/x = 2), we can derive closed- form expres- 
sions for V 2 : 
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C.4 The hyperbolic vacuum integral of the second 
kind 7^2 



C.4.1 Definition 

We define the hyperbolic vacuum integral of the second kind as 

3 r x 3 / 2 dx 
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and its domain is that of Tti(r, to). 
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C.4. 2 Properties 

The hyperbolic vacuum integral of the second kind is related to 7i\ (r, w) through 
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In the case of a flat (m = 1) perturbation, ?i 2 takes the form 
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which is the integral entering the linear growth factor in the Q m + £Ia = 1 universe. 
Hence, the linear growth factor function A(x) can be written as 

A(x) = 24/3( 3^ 2 2)1/2 ^(2~ 1/3 ^ 1) • (C18) 
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Appendix D 



Comparison to Bertschinger 



Similarity Solution 



In this appendix, we compare the Mach number as derived from the temperature jump 
across the accretion shock surrounding a collapsed structure (presented in chapter 3) 
to the similarity solution derived by Bertschinger [59] for the case of a single, spher- 
ically symmetric collapsed structure accreting matter in an otherwise homogeneous 
Q m = 1 universe. According to the Bertschinger solution, the accretion shock is posi- 
tioned at a constant fraction A s ~ 0.347 of the radius r ta of the matter shell turning 
around at a given cosmic time t, 



where m is the mass of the structure at the specific time t. The velocity of the infalling 
gas in the lab frame at the shock position is 
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The shock surface itself is propagating outwards, with a velocity v s = (8/9)A s r ta A l - 
Then, the absolute value of the gas velocity in the shock frame, v±, is 



Vi = \v g — v. 



■1.43- -X s 
9 



T -j (D.3) 



I.74 ( 1/3 (i + ^)i/2 m i/3 (D . 4) 



V 14 

If now the pre-shock material has an adiabatic sound speed c s i, the Mach number of 
the accretion shock is 

M(m,z,c sl ) = — 



1 + z Y' z m 



C S 1 

1.74 / 27GH \ 1/3 
I 14 j 

92 ( 151 ^ s " 1 )a+^) 1/2 



1/2^1/3 



> 1/3 



15 kms x \ u/2 / / c 



x — . (D.5) 

Note that the Bertschinger solution was derived in the limit A4 — ■> oo, and therefore 
it is only valid for M. 3> 1. 

We can compare Eq. (D.5) with the result we derived using the temperature jump 
across the surface of the shock. Combining Eq.s (3.3) (for the temperature jump in 
the high-AI limit) and (3.6) we get 

M « 93 ( 15 "'" )(l + 2 )'/ 2 

*(-)"". (D.8) 

in excellent agreement with Eq. (D.5). The small (~ 1%) deviation arises because, in 
obtaining M. from T2/T1, we have ignored any temperature structure within the col- 
lapsed object, and have instead assumed that the virial temperature of the structure 
is representative of the temperature right behind the shock, while the Bertschinger 
analysis calculates and takes into account the temperature structure inside the shock 
surface. 



1 if the time dependence of m is written out explicitly, the turnaround radius varies as rt a (i) oc t 8 / 9 
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